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ABSTRACT 

This  paper  presents  a  method  for  estimating  the  forces  and  moments  on  the  six 
axes  of  an  Unmanned  Underwater  Vehicle  (UUV)  operating  in  close  proximity  to  a  host 
vessel.  Forces  and  moments  encountered  by  the  vehicle  in  the  host  vessel  launchway,  near 
the  vessel  and  clear  of  the  launchway,  and  far  away  from  the  host  vessel,  are  modeled 
from  experimental  and  theoretical  sources.  A  quasi  -  steady  approach  to  the  vehicle 
launch  dynamics  allows  coupling  of  experimental  data,  potential  flow  theory,  and 
hydrodynamic  maneuvering  coefficients  to  develop  an  algorithm  capable  of  modeling 
vehicle  and  host  vessel  interactions  with  any  user  defined  trajectory. 

A  vehicle  motion  simulation  study  was  performed  using  the  Revised  Standard 
Equations  of  Motion  (EOM)  augmented  to  incorporate  the  interaction  force  model  for 
proximity  effects  of  the  UUV  to  a  host.  Application  of  the  trajectory  simulation  and  force 
and  moment  time  history  are  developed  for  a  underwater  vehicle  launching  from  a 
submerged  host  illustrating  the  hydrodynamic  interaction  of  the  vehicle  to  the  host.  The 
method  is  developed  to  account  for  host  forward  speed,  UUV  propulsor  effects,  and  open 
water  hydrodynamics.  The  robustness  of  the  method  is  demonstrated  for  several  launch 
trajectories. 

Thesis  Supervisor:  Dr.  Jerome  Milgram 
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Chapter  1    INTRODUCTION 

1.1        Motivation 

Underwater  vehicles  find  applications  in  an  increasing  number  of  diverse  missions 
and  roles.    Recent  uses  of  AUV's  and  ROV's  both  in  commercial  and  military 
applications  continue  to  outline  the  capabilities,  and  limitations  of  employing  unmanned 
underwater  vehicles  (UUV).  Integral  to  the  use  of  a  UUV,  whether  autonomous 
underwater  vehicle  (AUV)  or  remotely  operated  vehicle  (ROV),  is  the  process  used  for 
launch  and  recovery.  This  process  is  largely  driven  by  the  type  of  support  ship,  land  based 
launch  and  recovery  facility,  or  submerged  host  platform  Launch  and  recovery  from  a 
submerged  host  represents  an  opportunity  to  add  flexibility  to  the  vehicle  mission,  but 
adds  complexity  to  vehicle  design.  The  ability  of  a  submerged  host  vessel,  hereafter 
generically  referred  to  as  the  host,  to  deploy  a  ROV  or  AUV  to  a  mission  commencement 
location  in  latitude,  longitude  and  depth  offers  unique  advantages  in  mission  duration,  real 
time  interrogation,  and  covert  mission  accomplishment  over  other  launch  platforms.  The 
ability  to  properly  control  a  vehicle  as  it  transitions  from  within  the  submarine  launchway 
to  a  point  far  away  from  the  host  and  back  until  finally  rehoused  represents  a  complex 
dynamic  process  for  which  improved  models  are  necessary  to  investigate  the  forces, 
moments  and  trajectory  of  the  vehicle  during  the  course  of  its  mission  and  through  launch 
and  recovery.  The  design  of  effective  vehicle  control  is  the  goal  of  the  hydrodynamic 
modeling  of  the  launch  and  recovery  process. 
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1.2       Research  Objectives 

Launch  and  recovery  of  a  vehicle  from  a  submarine  is  expected  to  represent  the 
most  demanding  case  for  control  effectiveness.  The  effect  of  the  host  body  speed  and 
launchway  geometry  on  the  hydrodynamic  flow  characteristics  in  and  around  the 
launchway  have  been  documented  by  experiment  [11].  Additionally,  interaction  effects 
between  the  UUV  and  the  host  have  been  quantified  by  experiment  for  specific  cases  [11]. 
One  can  foresee  that  the  combined  effects  which  in  turn  are  coupled,  lead  to  a  complex 
hydrodynamic  modeling  problem 

The  first  objective  of  this  thesis  is  to  use  experimental  observations  from  a  scale 
UUV  within  and  in  the  vicinity  of  a  launchway  to  develop  a  model  to  describe  the  forces 
and  moments  acting  on  a  UUV  as  a  function  of  hydrodynamic  and  geometric  criteria 
relative  to  the  host.  Due  to  the  specific  nature  of  vehicle  and  host  geometry  on  the  forces 
experienced,  the  model  will  be  tailored  to  the  host  and  UUV.  However  this  specificity  is 
necessary  only  in  so  far  as  to  proceed  to  the  next  objective,  that  of  determining  a 
reasonable  simulation  to  determine  the  time  dependent  forces  and  moments  acting  upon 
the  vehicle. 

This  second  objective  is  to  incorporate  this  force  model  into  a  form  of  the  Revised 
Standard  Submarine  Equations  of  Motion  [5]  to  determine  the  forces  and  moments  acting 
upon  the  vehicle  given  a  specified  trajectory.  This  investigation  relies  upon  mating 
experimental  data  in  the  complex  region  within  the  launchway,  with  inviscid  flow  code 
results  external  to  the  launchway  and  in  proximity  to  the  host  "boundary",  and  accounting 
for  viscous  effects  of  cross  flow  drag  associated  with  the  vehicle  transitioning  through  a 
region  where  the  vehicle  axis  is  not  aligned  with  the  mean  free  stream  velocity.  The 
resulting  control  force  simulator  is  intended  to  be  adaptable  to  alternate  vehicle 
geometries. 


11 


Chapter  2    FORCE  MODEL 

2.1  Introduction 

The  forces  acting  upon  a  UUV  operating  in  the  open  ocean  differ  from  those 
experienced  when  operating  in  the  vicinity  of  a  boundary,  such  as  the  free  surface,  ocean 
bottom,  or  host  vessel.  For  the  purposes  of  this  thesis,  boundary  interaction  effects  are 
limited  to  a  plane  rigid  wall  approximation  of  a  host  vessel  except  in  the  internal  region  of 
the  launchway  of  the  host.  Interaction  forces  during  launch  and  recovery,  and  through  the 
transition  to  the  "infinite  fluid",  are  of  direct  interest  in  deterrnining  the  feasibility  of  the 
control  system  and/or  the  operating  envelope  for  launch  and  recovery.  Significant  work 
has  been  undertaken  in  the  past,  including  computational  modeling  of  these  forces  to 
assess  the  effects  of  high  speed  ejection  of  torpedoes  from  submarines  [2].  m  the  region 
of  interest  for  UUV  launch  and  recovery,  it  will  be  necessary  to  operate  within  a  different 
regime  of  slow  forward  speed  of  the  host  and  slow  relative  speed  of  the  UUV  to  the  host. 
This  combination  of  slow  host  and  UUV  relative  speed  indicates  the  plausibility  of 
modeling  the  UUV  launch  and  recovery  interaction  forces  as  quasi  -  steady,  the  UUV 
stationary  with  respect  to  the  host  and  slowly  changing  with  respect  to  its  variation  in 
geometric  relation  to  the  host.  This  approach  allows  use  of  experimental  results  to 
develop  a  model  for  the  forces  in  this  region.  Steady  flow  potential  based  codes  may  then 
be  used  to  model  the  interaction  effects  of  the  vehicle  near  the  host  to  resolve  areas  where 
experimental  data  is  incomplete. 

2.1.1  Background 

The  determination  of  the  hydrodynamic  forces  present  on  a  vehicle  transitioning 
from  a  submarine  launchway  is  complicated  by  the  number  of  parameters  involved  and  the 
complexity  of  the  fluid  -body  interactions.  For  a  baseline  case  with  no  vehicle  in  the 
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launchway  Sjoblom  and  Schwemin  [15]  showed  that  the  flow  in  a  typical  submarine 
launchway  is  unsteady  and  turbulent,  with  turbulence  intensities  within  the  cavity  ranging 
from  10  -  20%  of  free  stream  velocity.  Experimental  observations  document  average  flow 
magnitudes  0  -  10%  of  the  free  stream  at  the  aft  end  of  the  cavity  and  approximately  30% 
of  free  stream  at  the  midpoint  of  the  launchway  shutter.  They  document  a  large 
recirculation  region  in  mid  -  region  of  the  cavity  and  a  separated  region  approximately 
30%  of  the  cavity  length  forward.  As  reported,  the  turbulence  shows  a  flow  periodicity 
with  frequency  content  between  1  and  4  Hz.  External  to  the  launchway  and  turbulent 
boundary  layer  potential  flow  predictions  compare  well  with  measurements. 

In  order  to  assess  the  force  and  moment  characteristics  of  a  vehicle  in  launch  or 
recovery  transition,  both  the  region  within  the  launchway  and  external  to  it  require 
modeling.  The  model  must  be  sophisticated  enough  to  capture  three  dimensional 
geometry  effects,  body  -  to  -  body  interactions,  as  well  as  viscous  flow  effects.    Inviscid 
potential  flow  models  capture  the  proximity  interaction  effects;  however,  they  do  not 
include  viscous  drag  effects  including  cross  flow  drag  or  lift  forces  generated  by  control 
surfaces  and  vorticity  shed  into  the  wake  from  separation.  Numerical  potential  and 
viscous  flow  techniques  can  not  at  present  handle  all  these  effects  so  a  combination  of 
experimental  results,  numerical  techniques,  and  semi-empirical  non  linear  motion 
derivatives  are  suggested  for  investigating  the  launch  and  recovery  envelope. 

2.1.2  Experimental  Data 

Launchway  Zone 

Jensen  and  Milgram  [11]  describe  an  experiment  conducted  using  a  1/10  scale 
model  of  an  UUV  transitioning  through  a  submarine  launchway.  Figure  2-1  illustrates  the 
basic  configuration  and  geometry  associated  with  the  vehicle  and  the  host.  Principle 
dimensions  of  the  UUV  vehicle  are  listed  in  table  2-1. 
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UUV  Principle  Dimensions 

Length  (ft) 

10.5 

diameter  (ft) 

1.604 

Table  2- 1  Principle  Vehicle  Dimensions 


Figure  2- 1 
Data  was  obtained  for  the  forces  and  moments  at  various  points  of  advance  along  the 
launchway  centerline.  Reference  to  a  "guide  can"  structure  existing  in  the  torpedo  tube 
launchway  is  made  in  the  model  tests.  The  vehicle  was  advanced  from  a  position  farthest 
aft,  corresponding  to  the  aft  end  of  the  propulsor  duct  flush  with  the  forward  edge  of  the 
guide  can,  to  a  position  10  inches  forward  of  this  reference  (100  inches  full  scale).  This 
range  of  traverse  starts  with  the  vehicle  fully  within  the  launchway  and  extends  to  a 
position  with  the  aft  end  of  the  propulsor  duct  approximately  even  with  the  forward  most 
point  of  the  launchway.  Data  in  all  six  degrees  of  freedom  were  obtained  at  each  position 
along  the  launchway  for  different  onset  flow  conditions  corresponding  to  the  speed  of  the 
free  stream  past  the  host  allowing  the  effect  of  forward  speed  on  the  host  to  be  modeled. 
The  UUV  model  was  configured  with  a  ducted  propulsor  which  was  operated  during  the 
tests  at  various  RPM  for  each  advance  position  and  onset  flow  combination  to  simulate 
the  influence  of  propulsor  forces  and  race  effects  on  the  induced  flow. 
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Near  Boundary  Zone 

Once  the  vehicle  is  outside  of  the  launchway,  the  geometry  of  the  vehicle  and  the 
host  can  be  approximated  as  a  body  in  proximity  to  a  plane  rigid  boundary.  The  approach 
taken  here,  utilizes  results  from  a  Boundary  Integral  Element  Method  (BEEM)  panel  code 
to  determine  the  force  and  moment  effects  due  to  proximity  operations  external  to  the 
launchway.  This  method  uses  inviscid  potential  flow  theory  and  does  not  include  viscous 
effects.  A  description  of  the  method  and  implementation  is  contained  in  [14].  The  effect 
of  the  propulsor  on  the  induced  flow,  and  accompanying  modification  to  the  pressure 
integral  is  not  captured  by  the  BIEM  in  this  application.  An  attempt  to  include  these 
effects  is  made  through  use  of  experimental  data  obtained  for  the  UUV  at  various  onset 
flow  speeds,  propeller  RPM,  and  clearance  from  the  approximated  plane  boundary. 

2.2  Modeling  Results 

2.2.1  Coordinate  System 

The  force  and  moment  models  are  referred  to  a  right-hand  orthogonal  coordinate 
system  of  moving  axes,  fixed  in  the  UUV  with  the  origin  located  at  the  center  of  mass 
(CG)  of  the  vehicle.  The  positive  directions  of  the  axes  are  x  -  forward,  y  -  starboard,  and 
z  -  downward  in  accordance  with  Standard  Nomenclature  [16]. 

2.2.2  Launchway  Model 

The  experimental  data  lead  to  a  model  of  the  forces  and  moments  experienced  by 
the  vehicle  as  a  function  of  specific  dependent  variables.  The  parameters  varied  during  the 
testing  used  to  develop  the  model  were: 

(1)  onset  flow  speed 

(2)  advance  along  the  launchway  axis 

(3)  propulsor  rpm 
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(4)        thruster  flow  and  direction 
The  data  runs  with  the  thruster  running  were  not  used  in  this  analysis  in  order  to  limit  the 
force  and  moment  effects  to  those  induced  by  body  and  launchway  geometry  and 
propulsor  effects.  The  forces  and  moments  may  be  developed  in  accordance  with  the 
following  dimensional  form: 

F  =  f(d,U,p,x,y,z) 
M  =  f(d,U,L,p,x,y,z) 

The  following  non-dimensional  terms  are  defined  for  use  in  developing  a  model 

independent  of  scale. 

Advance 


CFY  = 


Length 
F 


\pVof2Af 


CMZ  = 


IpUjAjL 


___,     (  PROPULSIVE  FORCE  ) 
Fib  = 

~2pVof2Af 


where:    CF  =  Force  Coefficient 

CM  =  Moment  Coefficient 
PTF  =  Propulsive  Thrust  Factor 
Uof  =  Onset  flow  speed 
Ay  =  Frontal  Area 
L  =  Vehicle  Length 


For  each  data  set  the  forces  and  moments  recorded  were  converted  to  standard  body 
coordinates  with  forces  and  moments  limited  to  the  horizontal  plane  being  modeled  in  this 
case.  A  propulsion  model  was  required  which  would  allow  robust  modeling  of  propeller 
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effects.  The  non-dimensional  Propulsive  Thrust  Factor  (PTF)  was  used,  with  PTF 
determined  for  each  data  set  by  measuring  the  mean  recorded  thrust  at  zero  RPM  for  each 
onset  flow  condition  and  advance  to  length  position  and  subtracting  this  value  from  the 
axial  force.  The  zero  PTF  data  corresponds  to  the  axial  drag  force  experienced  in  the 
launchway  due  to  the  body  position  alone  and  was  modeled  as  well.  The  procedure 
selected  to  arrive  at  a  model  recognized  that  significant  variance  existed  in  the 
experimental  data.  Figures  2.2  through  2. 14  reflect  the  non  dimensional  data  with 
associated  best  least  squares  linear  regression  fit. 
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Figure  2.2  Drag  Force 
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CFY(PTF)  for  A/L  0.68 


Figure  2.3  Lateral  Force  for  A/L  =  0.68 

CFY(FTF)  for  A/L  0.84 
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Figure  2.4  Lateral  Force  A/L  =  0.84 
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CFY(FTF)  for  A/L  1.0 


t 


Figure  2.5  Lateral  Force  A/L  =  1.0 

CFY(PTF)  for  A/L  1.16 


Figure  2.6  Lateral  Force  A/L  =  1. 16 
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CFY(FTF)  for  A/L  1.32 
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Figure  2.7  Lateral  Force  A/L  =1.32 
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Figure  2.8  Lateral  Force  A/L  =1.48 
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CMZ(PTF)  for  A/L  0.68 
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Figure  2.9  Lateral  Moment  A/L  =  0.68 
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Figure  2. 10  Lateral  Moment  A/L  =  0.84 
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CMZ(FTF)  for  A/L  1.0 


Figure  2. 1 1  Lateral  Moment  A/L  =1.0 

CMZ(PTF)  for  A/L  1.16 


Figure  2. 12  Lateral  Moment  A/L  =  1. 16 
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CMZ(PTF)  for  A/L  1.32 


Figure  2. 13  Lateral  Moment  A/L  =1.32 

CMZ(PTF)  for  A/L  1.48 


Figure  2. 14  Lateral  Moment  A/L  =  1.48 
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A  One  Way  Analysis  of  Variance  (ANOVA)  was  conducted  on  each  data  set  of  force  and 
moment  as  a  function  of  A/L  and  PTF.  Table  2-2  summarizes  the  results. 


A/L 

Data 

F  Ratio 

Prob>F 

r2 

0.68 

CFY 

4.0022 

0.0627 

0.200087 

CMZ 

539.6004 

<0.0001 

0.971202 

0.84 

CFY 

2.6603 

0.1339 

0.210132 

CMZ 

2.2577 

0.1639 

0.184186 

1.0  f 

CFY 

276.3695 

<0.0001 

0.802538 

CMZ 

44.7450 

<0.0001 

0.396869 

1.16 

CFY 

88.9794 

<0.0001 

0.847589 

CMZ 

167.6809 

<0.0001 

0.912892 

1.32 

CFY 

166.9032 

<0.0001 

0.830764 

CMZ 

116.0631 

<0.0001 

0.773429 

1.47 

CFY 

31.1761 

<0.0002 

0.757141 

CMZ 

857.6477 

<0.0001 

0.988475 

f  (2)  outlying  points  omitted 
Table  2-2  One  Way  ANOVA  Results 
In  Table  2-2  ,  the  F-ratio  is  a  statistical  measure  indicating  whether  the  model  is  a 

better  fit  than  the  response  mean,  as  determined  by 

(V\,     Degrees  of  freedom) 

[/ \yt     Degrees  of  freedom) 


F  -  ratio- 


model 


The  r^  category  represents  the  amount  of  total  variation  accounted  for  by  the  model  with 
the  remaining  variation  attributed  to  random  error  and  is  determined  as 


model 


CtotaJ 


To  conclude  from  Table  2-2  that  the  force  and  moment  are  dependent  upon  PTF  the  F- 
ratio  must  exceed  an  observed  significance  probability  (P- Value)  of  0.05.  To  determine 
the  effectiveness  of  the  model,  r^  is  used  in  conjunction  with  the  F-ratio.  The  case  for 
A/L  of  1.0  utilized  the  r^    parameter  to  omit  two  outlying  points  from  the  model.  The 
results  indicate  that  for  the  majority  of  sample  conditions  the  effect  of  PTF  at  each  A/L  is 
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present,  and  that  a  linear  model  is  justified.  Higher  order  models  were  not  justified  In 
order  to  arrive  at  a  smooth  model  for  force  and  moment  as  a  function  of  PTF  and  A/L  the 
cCfy{ptf)Jcptfa  f  L   is  used  to  fit  a  3rd  order  cubic  function  and  the  intercept  is  adjusted  to 

make  a  best  fit.  Figure  2-15  and  2-16  reflect  the  best  fit  for  force  and  moment  within  the 
launchway.  It  can  be  seen  that  a  positive  outboard  force  exists  for  a  minimum  A/L  of  0.64 
and  remains  positive  throughout  the  vehicle  transition  along  the  launchway  centerline  with 
a  maximum  value  at  an  A/L  of  1. 15  which  corresponds  to  approximately  one  quarter  of 
the  stern  shadowed  by  the  launchway.  The  lateral  force  then  diminishes  to  50%  of  this 
force  at  the  maximum  A/L  of  1.48.  Lateral  moment  effects  are  quite  different,  showing  an 
inboard  moment  at  minimum  A/L  followed  by  a  zero  crossing  at  A/L  =  0.9.  From  this 
position  the  moment  shows  a  continual  increase  to  a  maximum  with  the  vehicle  fully 
emerged  from  the  host,  corresponding  in  magnitude  but  opposite  sign  its  value  at 
minimum  A/L. 
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5CFY(PTF)/5PTF  VS  A/L 


Figure  2.15  Fit  of  xfy \ptf)j 'cPtf a  f  L 

3CMZ(PTF)/ePTF  VS  A/L 


Figure  2. 16  Fit  of  cCmz(Ptf)/cPtfa  f  L 
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The  resulting  force  and  moment  model  in  the  launchway  are  expressed  as 

CFY  =  fi(sCFY(PTF)/d-TFA/  L)A  /  L  +  /.(CFY^^) 

CMZ  =  f3(cCMZ(PTF)/cPTFA ,  l)a  I  L  +  fx{CM2 Mercepl) 
where  fx  represents  a  linear  fit 

f^  represents  a  third  order  polynomial  fit 
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Figure  2. 17  Sample  of  Lateral  Force  Model  Fit  A/L  =  1.0 
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CMZflTF)  for  A/L  1.0 
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Figure  2. 18  Sample  of  Lateral  Moment  Model  Fit  A/L  =1.0 

2.2.2  Host  Proximity 

With  the  vehicle  near  the  boundary  but  not  in  the  launchway,  insufficient  data 
existed  in  [1 1]  to  predict  the  effects  of  relative  pitch  and  clearance  from  the  boundary. 
To  develop  this  model,  an  inviscid  potential  flow  panel  code  was  resorted  to,  in  order  to 
determine  the  near  host  interaction  effects.  Difficulties  in  the  numerical  hydrodynamics 
were  encountered  in  modeling  due  to  the  blunt  forward  end  of  the  UUV  geometry.  Initial 
paneled  models  attempted  to  reflect  the  bow  accurately,  however  the  large  change  in 
sectional  area  vs.  length  resulted  in  a  large  negative  pressure  at  the  bow  leading  to 
excessive  axial  suction  force  instead  of  the  zero  drag  condition  expected.  Modification  of 
the  model  to  blend  the  change  in  bow  cross  section  over  increased  length  provided  a 
model  which  corresponded  well  with  strip  theory  and  experimental  data  to  a  similar 
vehicle  [10].  The  resulting  panel  geometry  is  illustrated  in  Figure  2. 19. 
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Figure  2. 19  UUV  Panel  Geometry 
The  proximity  model  is  obtained  in  two  parts.  First,  investigating  the  effect  of  relative 
proximity  to  the  boundary  both  with  and  without  relative  pitch  including  far  field 
determination  of  the  Munk  moment  using  potential  flow.  And  second,  investigating  the 
effect  of  PTF  on  force  and  moment  from  experimental  data.  Figure  2.20  shows  the 
suction  force  and  moment  as  a  function  of  gap/diameter,  where  gap  refers  to  the  distance 
from  the  vehicle  molded  shell  to  the  boundary.  The  effect  of  pitch  and  gap/diameter  on 
the  force  and  moment  is  illustrated  for  selected  cases  in  Figures  2.21  through  2.23.  The 
non  linear  least  square  fit  is  superimposed  in  Figures  2.21  and  Figure  2.22  showing  good 
agreement  over  the  range  of  combinations  of  pitch  and  gap/diameter.  Figure  2.23 
represents  the  potential  flow  Munk  moment  determined  with  the  vehicle  in  an  unbounded 
fluid  domain. 
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CFY  and  CMZ  At  Zero  Fitch 
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Figure  2.20  Boundary  Interaction  At  Zero  Pitch 


CFY  vs  gap/dia  for  Pitch  From  0°  to  20° 
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Figure  2.21  Lateral  Force  Near  a  Boundary 
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CMZ  vs  gap/diameter  for  Pitch  of  1° 
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Figure  2.22  Lateral  Moment  Near  a  Boundary 
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Figure  2.23  Open  Water  Munk  Moment 
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The  effect  of  PTF  on  the  force  and  moment  is  captured  by  data  obtained  in  [1 1].  Figures 
2.24  and  2.25  represent  the  lateral  force  and  moment  effects  due  to  a  change  in  propulsive 
thrust  factor.  Significant  scatter  of  experimental  data  was  observed  and  several  outlier 
points  have  been  omitted  in  both  curve  fits,  using  the  method  described  previously,  to 
arrive  at  a  model  which  makes  physical  sense  and  tends  to  the  proper  limit. 
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Figure  2.24  PTF  Effect  on  Lateral  Force 
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CMZ(PTF)  Near  Wall 
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Figure  2.25  PTF  Effect  on  Lateral  Moment 
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Chapter  3    SIMULATION  MODEL 

3.1  Coordinate  Systems 

The  forces  acting  on  an  underwater  vehicle  are  normally  referred  to  a  non-inertial 
coordinate  system  fixed  in  the  body,  however,  motion  of  the  vehicle  is  most  easily 
described  in  an  inertial  reference  frame.  For  the  purposes  of  this  analysis,  the  inertial 
reference  frame  is  translating  at  the  constant  velocity  of  the  host,  assumed  to  be  along  the 
positive  X  axis,  and  therefore  represents  a  host  fixed  reference.  The  Z  axis  is  directed 
positive  downwards  and  the  positive  Y  axis  is  to  starboard  in  this  frame  of  reference.  To 
obtain  the  forces  and  moments  in  the  non-inertial  reference  it  is  necessary  to  transform  the 
vehicle  motion  to  this  frame.  The  vector  transformation  between  the  reference  frames  is 
uniquely  determined  by  three  angles:  $  rotation  about  the  X  axis  (roll),  0  rotation  about 
the  Y  axis  (pitch),  vjy  rotation  about  the  Z  axis  (yaw),  which  describe  the  angular 
orientation  between  the  two  frames  and  are  termed  the  Euler  angles.  The  order  in  which 
the  transformation  about  these  angular  relations  is  observed  is  significant,  and  if  taken  in 
different  sequence  can  produce  a  different  vector  result.    In  accordance  with  standard 
convention  the  transformation  is  accomplished  by  ordering  the  rotations  as  first  yaw  ( \j/), 
then  pitch  (0),  and  finally  roll(<j)).  The  transformation  matrix  describing  the  vehicle  motion 
in  translation  in  the  non-inertial  frame,  given  the  inertial  frame  motion,  is  expressed  by: 


Tii*,0,r) 


cos(#)cos(^)  cos(#)sin(^)  -sin(#) 

-  sin(  y/)  cos(  </>)  +  sin(  <p)  sin(  6)  cos(  y/)      cos(  </>)  cos(  y/)  +  sin(  tj>)  sin(  6)  sin(  y/)      sin(  <p)  cos(  0) 
sin(  </>)  sin(  y/)  +  cos(  <j>)  cos(  y/)  sin(  0)      -  sin(  tp)  cos(  y/)  +  cos(  <j>)  sin(  6)  sin(  y/)    cos(  (p)  cos(  0) 


This  transformation  accounts  for  the  rigid  body  rotation  from  host  fixed  reference  to  body 
fixed  reference.  The  transformation  of  angular  rate  quantities  follows  the  same  Euler  angle 
transformation  convention  and  order  of  rotations. 
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The  angular  rate  transformation  is  then: 


T2W,07V,)  = 


1  0  -sin(<?) 

0     cos(^)      sin(^)cos(^) 
0    -sin(^)    cos(#)cos($) 


The  total  transformation  matrix  may  be  represented  as: 


T{<f>,d,y)  = 


Tx^,e,yf)  0 

0  T2{<j>^¥) 


With  <ft0  representing  a  vector  in  the  inertial  frame  consisting  of  scalar  quantities  of  u 
(surge),  v  (sway),  and  w  (heave)  velocities  ,  and  p  (roll  rate),  q  (pitch  rate)  and  r  (yaw 
rate),  or  their  respective  time  derivatives;  the  previously  defined  Euler  angles  (|),  0,  and  \\j 
are  used  to  determine  the  vector  <K,  representing  the  same  vector  SR0,  but  in  the  body 
reference  frame  changed  only  in  orientation.  The  relation 

is  used  to  determine  the  body  relative  velocity  and  acceleration  in  the  body  reference 
frame  translating  with  the  velocity  of  the  host  in  the  host  axial  direction.  Finally,  the  body 
velocities  are  corrected  for  this  constant  translation  of  the  body  fixed  reference  frame  to 
arrive  at  the  total  body  relative  velocities  independent  of  the  host. 

3.2  Equations  of  Motion  and  Implementation 

The  equations  of  motion  of  a  rigid  body  may  be  expressed  for  arbitrary  origin  in 
matrix  form  as  follows: 
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This  expresses  the  Revised  Standard  Equations  of  Motion  [5]  with  the  acceleration  terms 
isolated  to  the  left  hand  side.  The  right  hand  side  consists  of  the  sum  of  the  gyroscopic 
and  centripetal  acceleration  effects  as  well  as  forces  arising  from  the  fluid  hydrodynamic 
effects  on  the  vehicle.  These  last  include  forces  and  moments  accounted  for  by 
maneuvering  coefficients  including  control  forces,  hydrostatics,  crossflow  and  vortex 
effects.  The  hydrodynamic  forces  are  calculated  in  accordance  with  the  standard 
equations  of  motion.  These  are  described  for  the  hydrodynamic  terms  as: 


X  =  Xqqq2  +  Xrrr2  +  Xrprp  +  Xvrvr  +  Xwq  wq  +  Xmv2  +  Xwvw2  +  Xwuw  +  Xquq  +  Xuuu2 
Y  =  yp\P\P\p\  +  YpqPq  +  Yrur  +  Ypup  +  Yvuv  +  Yv]^v\uv\  +  Y„v2 +Ywpwp  +  Yvluvy\uv\  +  YMr\r\  +  Y,uu  +  Yv{rv{v 
Z  =  Zquq  +  ZwUW  +  Z^u\w\  +  Z^vp  +  Z.uu  +  Zq^q\q\  +  Zw^w\w\  +  Zrprp  +  Zw  w 
K  =  Kpup  +  Krur  +  Kqrqr  +  K^wp  +  K^w  +  Kvuv  +  KpMp\p\  +  Kmu2 
M  =  Mquq+Mq^q\q\  +  M^w\w\  +  Mvpvp+Mrpr  +  Mmu2  +Mwuw  +  Mwww-<Jv2  +w 

4J~+ 


I    2    ,       2 
\V     +  W 


N  =  N    pq  +  N  up  +  N    wp  +  N.u*  +  Nrur  +  Nvuv  +  Nv{vr{r 


+  NMu\v\  +  Nwv2  +  NMv\v\  +  N&M 


The  body  terms,  including  coriolis,  centripetal,  and  hydrostatic  effects  are  calculated 
separately.  Neglecting  cross  products  of  inertia  these  terms  may  be  expressed  in  the  form: 

Xbody  =-n\wq  -vr-xG(q2+r2)+yGpq  +  zG(pr))-(W-B)sin(ff) 

Ybody  =-m{ur-wp  +  xGpq  +  zGqr)  +  (W  -  B)cos(ff)sin(4>) 

Zbody  =-nAyp-^]  +  xGpr-zG{p2  +  g2)  +  (fl/-5)cos(6>)cos(<z>) 

Kbody  =(Iyy-I2z)qr-m{yG(vp-wq)  +  zG(ur-\ip))  +  (yGW-yBB)cos(0)cos(^)-(zGfV-zBB)cos(0)sin(<f>) 

Mbody  =(Izz-Ixx)pr-m(xG(\T-uq)-zG(wq-w))-(xGW-xBB)cos(0)cos(<fi)-(zGW-zBB)sm(ff) 

Nbody  =(Ixx  -Iyy)pq-m{xG(ur-Hp)-yG(wq  -vr))  +  (xGW-xBB)cos(0)sm(<t>)  +  (yGtV-yBB)sin(0) 
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The  hydrodynamic  maneuvering  coefficients  may  be  obtained  by  experiment,  or  predicted 
by  numerous  methods  [1],  [9],  [12].  These  relations  do  not  include  the  propulsive  thrust 
model  forces  or  control  forces  in  order  to  prevent  redundancy  of  these  terms  with  the 
interaction  force  model. 

The  right  hand  side  forces  on  the  vehicle  may  be  expressed  as  a  partition  of  the 
total  right  hand  side  in  the  following  fashion: 

Ftou.  =  Fh  +  Fi  +  Fc 
where  Ftota,  may  be  equated  to  the  acceleration  and  mass  product  terms  ,  FH  represents  the 

forces  and  moments  acting  on  the  body  due  to  open  water  hydrodynamic,  centripetal, 
gyroscopic,  and  hydrostatic  effects.  The  term  Fj  accounts  for  the  interaction  effects  due  to 
the  operation  in  a  restricted  fluid  e.g.  in  the  launchway  and  host  vicinity.    Fc  accounts  for 
the  control  force  and  moment  required  to  effect  the  selected  trajectory.    Fc  may  be  solved 
for  directly: 

Fc  =  Ft0Ul  -  (FH  +  Fj) 
This  relation  is  adapted  for  implementation  in  a  six  degree  of  freedom  motion  simulation 
program  for  which  a  general  description  of  the  method  employed  follows. 

The  desired  trajectory  of  the  vehicle  must  be  known  a  -  priori  and  specified  as 
input  to  the  algorithm  It  is  important  in  the  launch  and  recovery  of  the  vehicle  to  know 
the  relative  position  ,  velocity,  and  acceleration  of  the  vehicle  with  respect  to  the  moving 
host,  therefore  the  trajectory  is  specified  in  a  host  fixed  inertial  frame.  The  vehicle 
accelerations  and  velocities  are  determined  by  a  coordinate  transformation  from  the 
moving  inertial  frame  to  the  body  fixed  frame  with  the  effect  of  forward  speed  of  the  host 
accounted  for  in  the  transformation.  The  total  force  at  each  time  step  is  determined  from 
the  vehicle  accelerations  by  solving  the  product  of  the  mass  matrix  and  body  relative 
acceleration  vector  directly.  A  word  of  explanation  is  warranted  with  respect  to  the  use  of 
a  motion  simulation  program  For  the  aspect  of  this  thesis,  with  the  specification  of  a 
trajectory  known  explicitly  over  time  along  with  acceleration,  and  velocity,  there  is  no 
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need  to  invert  the  mass  matrix  and  integrate  over  each  time  step.  The  forces  and  moments 
may  be  computed  directly  at  each  time  step  allowing  the  algorithm  to  run  very  fast  without 
growing  errors  due  to  a  integration  scheme.  The  hydrodynamic  forces  and  moments 
including  the  body  terms  are  calculated  at  each  time  step  without  consideration  to  the  host 
presence.  A  separate  subroutine  implements  the  interaction  force  model  previously 
developed,  where  forces  and  moments  when  the  vehicle  is  located  within  the  launchway 
are  determined  completely  from  the  launchway  model  and  include  hydrodynamic  and 
interaction  effects.  Once  the  vehicle  has  transitioned  to  a  position  outside  of  the 
launchway  the  interaction  effects  alone  are  determined.  This  condition  reflects  only  a 
portion  of  the  forces  acting  upon  the  vehicle,  therefore  to  arrive  at  the  total  right  hand  side 
force  vector  when  the  vehicle  is  outside  of  the  launchway,  the  hydrodynamic  derivative 
based  forces  are  added  to  the  interaction  effects. 

In  order  to  achieve  a  smooth  transition  between  the  experimental  data  representing 
the  total  force  on  the  vehicle  to  the  force  calculated  by  the  hydrodynamic  derivative  based 
force  vector,  it  is  necessary  to  employ  a  smoothing  procedure  to  blend  the  differences 
between  the  two.  In  order  to  accomplish  this,  a  zone  structure  is  defined  with  the 
following  description: 

Zone  1     Inside  the  host  launchway  F(x,y,z)  =  Hydrodynamice  +  Interactione 
Zone  2     Within  the  transition  region  F(x,y,z)  =  Hydrodynamic^  +  Interaction^ 
Zone  3     Outside  of  the  transition  region  F(x,y,z)  =  Hydrodynamic^  +  Interactione 
Where  F(x,y,z)  is  the  force  calculated  by  the  force  model  and  subscripts 
e>  ^  and  d  refer  to  experiment,  interpolated,  and  coefficient  based 
semi- empirical  data  respectively. 

Zone  2  is  where  the  smoothing  routine  operates.  This  zone  is  defined  by  a  hemisphere 
with  the  center  located  at  the  launchway  exit  defined  by  xc,yc,zc  with  the  diameter  of  this 
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4)  A  geometry  limitation  file  containing  the  position  information  with  regard 
to  launchway  geometry,  host  geometry  and  host  speed. 

5)  A  trajectory  file  which  includes  vehicle  acceleration  and  velocity  scaling 
parameters  and  trajectory  modification  commands  versus  time. 

A  linear  interpolation  scheme  is  used  to  determine  the  body  position,  velocity  and 
acceleration  changes  over  time  from  the  commanded  trajectory  inputs.  A  partial  list  of  the 
subroutines  implementing  the  trajectory,  force  model  and  coordinate  transformation  is 
included  as  Appendix  A. 
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Chapter  4    SIMULATION  RESULTS 

4.1  Trajectory  Runs 

Four  different  trajectories  were  selected  for  evaluation  of  the  simulation.  A  copy 
of  the  input  files  is  included  in  Appendix  B.  The  trajectories  were  selected  to  test  the 
blending  between  the  force  and  moment  models  between  Zone  1  and  Zone  3.  The 
trajectories  selected  were  not  intended  to  reflect  in  all  cases  the  trajectory  that  would  be 
optimal  or  would  represent  the  real  vehicle  desired  course  but  would  illustrate  the  model. 

4.2  Straight  Launch 

The  first  trajectory  simulated  was  a  straight  ahead  motion  of  the  vehicle.  Using  an 

exponential  time  function  acceleration  described  by: 

x(t)  =  2(l-e-*-°*,i) 

(-0.04t2) 

xBody(t)  =  (4*0.04teK  ') 

The  vehicle  is  directed  out  of  the  launchway  and  away  from  the  host  at  a  constant  heading 
angle.  The  trajectory  in  host  reference  is  illustrated  in  Figure  4.1.  Figure  4.2  through  4.7 
show  the  resulting  accelerations,  velocities,  forces  and  moments  for  this  maneuver  in  the 
vehicle  reference  frame. 
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Figure  4. 1  Vehicle  Trajectory 
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Figure  4.2  Velocity  Profile 
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Traj_l  Vehicle  Relative  Acceleration 
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Figure  4.3  Acceleration  Profile 
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Figure  4.4  Hydrodynamic  Forces  in  Open  Water 
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Trajl  Interaction  and  Hydrodynamic  Forces 
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Figure  4.5  Hydrodynamic  Forces  Including  Vehicle  Interactions 
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Figure  4.6  Total  Force 
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Traj_l  Control  Force 
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Figure  4.7  Control  Forces 
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Although  unrealistic  as  a  desired  maneuver  due  to  the  large  sustained  crossflow  drag 
resulting  from  the  constrained  yaw  angle,  this  trajectory  is  illustrative  of  the  interaction 
effects  and  transition  from  a  region  where  interaction  effects  play  a  role  to  the  region  in 
the  far  field  where  the  interaction  effects  are  small.  As  such  it  is  useful  as  a  baseline 
condition.  The  location  of  the  vehicle  relative  to  the  zones  previously  defined  is  indicated 
in  the  figures,  where  pertinent,  by  a  zone  indicator  line  which  takes  a  value  of  1,  2,  or  3 
respective  to  the  appropriate  zone. 

The  interaction  effects  of  the  launchway  and  host  vehicle  are  apparent  in  the  figures.  As 
the  vehicle  nose  passes  the  launchway  exit  (occurring  very  early  in  the  trajectory  at  ~  4 
sec)  the  lateral  force  can  be  seen  to  decrease  as  the  nose  proceeds  into  the  transition  zone. 
The  effect  on  the  lateral  moment  indicates  an  almost  impulsive  effect  turning  the  moment 
from  a  negative  value  when  the  vehicle  is  fully  within  the  launchway  to  a  large  positive 
moment  as  the  vehicle  fully  exits.  The  magnitude  of  this  moment  depends  upon  the 
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interaction  due  to  the  boundary  as  well  as  the  effect  of  PTF  on  the  accelerated  flow 
between  the  vehicle  and  the  interacting  boundary.  As  the  boundary  recedes,  the  moment 
decreases  to  the  open  water  Munk  moment  due  to  the  imposed  angle  of  attack  of  the 
vehicle  in  this  trajectory  simulation. 

4.3  Launch  with  Yaw  Maneuver 

The  second  trajectory  for  simulation  yaws  the  vehicle  to  a  zero  relative  heading 
angle  after  the  vehicle  has  cleared  the  launchway.  This  represents  a  more  realistic 
maneuver  intended  to  investigate  the  effect  of  maneuvering  to  remove  crossflow.  The 
vehicle  is  yawed  to  a  zero  heading  angle  at  time  {t  =  25  sec}  and  proceeds  along  a  path 
parallel  to  the  host  for  the  rest  of  the  simulation.  This  trajectory  is  illustrated  in  Figure 
4.8.  Figures  4.9  through  4. 14  illustrate  the  velocities,  forces  and  moments. 


Traj_2  Position  History 


300 


200 


100 


V[deg] 
y  position 
x  position 


300 


200 


■  100 


20 


40 


60 


80 


— '0 
100 


Time  [s 


Figure  4.8  Vehicle  Trajectory 
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Traj_2  Body  Relative  Velocities 
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Figure  4.9  Velocity  Profile 
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Figure  4. 10  Acceleration  Profile 
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Traj_2  Open  Water  Hydrodynamic  Force 
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Figure  4. 1 1  Hydrodynamic  Forces  in  Open  Water 
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Figure  4. 12  Hydrodynamic  Forces  Including  Vehicle  Interactions 
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Figure  4. 13  Total  Force 
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Figure  4. 14  Control  Force 
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The  interaction  effects  near  the  host  are  apparent  in  the  lateral  force  and  moment  as  they 
were  in  the  straight  launch  trajectory.  In  addition,  acceleration  effects  during  the 
maneuver  can  be  seen  to  affect  the  required  control  forces,  specifically  in  this  instance  the 
contribution  of  N-  to  the  lateral  moment. 


4.4  Launch  with  Yaw  and  Deceleration  Maneuver 

The  third  trajectory  selected  for  simulation  illustrates  a  maneuver  intended  to 
reflect  a  realistic  vehicle  maneuver.  In  this  simulation  the  vehicle  transitions  out  of  the 
launchway  and  subsequently  decelerates  to  the  desired  speed  once  clear  of  the  host. 
Figure  4. 15  shows  the  trajectory  history.  Figures  4. 16  through  4.21  illustrate  the 
velocities,  forces  and  moments. 
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Figure  4.15  Vehicle  Trajectory 
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Figure  4. 16  Velocity  Profile 
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Figure  4. 17  Acceleration  Profile 
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Traj_3  Open  Water  Hydrodynamic  Force 
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Figure  4. 18  Hydrodynamic  Forces  in  Open  Water 
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Figure  4. 19  Hydrodynamic  Forces  Including  Vehicle  Interactions 
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Figure  4.20  Total  Force 
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Again  as  in  the  previous  trajectory  the  significance  of  the  acceleration  effects  are  evident 
in  the  axial  and  lateral  body  forces  during  deceleration.  The  lateral  moment  peak  occurs 
just  after  the  vehicle  has  cleared  the  launchway  and  experiences  the  maximum  interaction 
effect  of  the  planar  boundary  in  proximity  to  the  vehicle  in  this  configuration. 

4.5  Straight  Launch  At  Zero  Host  Speed 

The  fourth  trajectory  selected  for  simulation  illustrates  a  launch  from  a  stationary 
host.  In  this  simulation  the  vehicle  transitions  out  of  the  launchway  maintaining  a  straight 
course  and  speed  without  onset  flow.  Figure  4.22  shows  the  trajectory  history.  Figures 
4.23  through  4.28  illustrate  the  velocities,  forces  and  moments. 
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Figure  4.22  Vehicle  Trajectory 
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Traj_4  Body  Relative  Velocities 
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Figure  4.23  Velocity  Profile 
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Figure  4.24  Acceleration  Profile 
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Traj_4  Open  Water  Hydrodynamic  Force 
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Figure  4.25  Hydrodynamic  Forces  in  Open  Water 
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Figure  4.26  Hydrodynamic  Forces  Including  Vehicle  Interactions 
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Figure  4.28  Control  Force 
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It  is  of  interest  to  note  the  differences  between  this  launch  and  the  previous  launches,  with 
forward  host  speed.  In  this  launch  the  interaction  effects  are  clearly  seen;  the  lateral  force 
changing  sign  from  an  outward  force  when  the  vehicle  is  mostly  within  the  launchway  to  a 
inward  suction  force  once  clear.  The  lateral  moment  also  changes  sign  from  a  bow  in 
moment  when  the  vehicle  is  mostly  within  the  launchway  to  an  outward  moment  due  to 
boundary  interaction  and  suction  region  acting  on  the  aft  portion  of  the  vehicle.  The  axial 
force  reflects  the  acceleration  to  constant  velocity  with  corresponding  constant  axial 
hydrodynamic  drag,  while  the  lateral  moment  approaches  asymptotically  to  a  value 
dependent  upon  PTF.  The  effect  of  forward  host  speed  has  a  marked  effect  on  the  force 
and  moment  on  the  vehicle,  from  commencement  of  launch  with  the  vehicle  fully  within 
the  launchway,  to  its  transition  to  open  water.  With  a  4  feet  per  second  onset  flow 
velocity  due  to  the  host  forward  speed,  the  lateral  force  does  not  achieve  sign  reversal, 
and  an  outward  force  dominates  throughout  the  maneuver  sequence.  The  effect  on 
moment  is  pronounced  as  well,  reflecting  a  twenty  fold  increase  in  peak  moment  over  the 
zero  forward  speed  case. 

The  small  perturbation  of  the  lateral  moment  observed  as  the  vehicle  becomes  fully 
emerged  from  the  launchway  is  qualitatively  similar  to  the  moment  profile  results 
described  in  [2]  for  a  high  speed  torpedo  ejection  from  a  submarine  launchway.  This 
comparison  is  not  intended  to  validate  the  results  here  but  is  offered  for  completeness. 
Due  to  the  disimilar  nature  of  the  ejection  dynamics  compared  to  vehicle  swimout,  a 
quantitative  comparison  was  not  conducted. 
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Chapter  5    CONCLUSION 

5.1  Summary 

From  the  comparison  of  trajectory  runs  it  is  clear  that  the  force  and  moment  model 
implemented  into  a  six  degree  of  freedom  motion  algorithm  can  solve  for  control  forces 
versus  time  for  any  desired  trajectory.  The  mathematical  model  successfully  makes  non 
dimensional  and  integrates  experimental  data  with  potential  flow  panel  code  results,  and 
with  a  parametric  based  non  linear  mathematical  model,  to  span  the  trajectory  regime  from 
launchway  to  open  water.  Host  vessel  speed  effects  are  incorporated.  The  smoothing 
technique  accounts  for  the  inevitable  mismatch  between  the  experimentally  derived 
interaction  model  to  the  parametric,  coefficient  based,  hydrodynamic  forces  in  open  water 
and  allows  for  geometrically  similar  models  to  be  tested  using  the  existing  experimental 
results. 

Based  upon  trajectories  with  changing  accelerations,  acceleration  effects  add 
significant  information  to  the  force  simulation.  Vehicle  motion  gain  evaluation  is  enabled 
through  the  trajectory  input  file. 

The  magnitudes  of  the  required  control  forces  are  broken  down  in  such  a  way  so 
as  to  enable  the  contributing  interaction  forces,  open  water  hydrodynamic  forces,  or  body 
mass  forces,  to  be  compared  and  their  dependence  upon  trajectory  evaluated. 

5.2  Recommendations  for  Future  Work 

The  control  force  output  from  the  model  represents,  in  a  sense,  the  magnitude  of 
control  actuated  force  to  maintain  the  selected  trajectory.  A  shortcoming  in  the  model 
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exists  in  that  it  ignores  interaction  effects  between  the  control  surface  or  thruster,  the 
vehicle,  and  proximity  boundary.  The  fidelity  of  the  simulation  would  be  improved  by 
incorporating  control  effects  on  the  uncontrolled  vehicle  interaction  force  model. 

For  a  vehicle  with  a  different  propulsor,  a  means  of  effectively  adapting  the  thrust 
model  and  resulting  thrust,  moment  and  torque  effects  would  be  useful. 

The  basis  for  a  maneuvering  workbench  exists  which  can  be  used  to  assist  in 
control  appendage  sizing,  control  system  design,  and  free  swimming  motion  simulation  in 
one  simulation  environment. 

The  accuracy  of  the  simulation  and  its  representation  of  the  real  world  can  only  be 
determined  by  comparing  its  output  to  full  scale  tests.  These  would  indicate  those  areas 
captured  by  the  model  as  well  as  areas  requiring  improvement.  Additional  scale  testing 
should  be  conducted  to  provide  a  more  robust  mathematical  model,  and  should 
incorporate  thruster,  propeller,  and  launchway  off-  axis  effects. 
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Appendix  A.  Source  Code  Listing  for  Trajectory  Simulation  and  Force  Model 

This  appendix  includes  source  code  for  the  implementation  of  the  trajectory,  and  force 
model.  Also  included  is  a  portion  of  the  Equation  of  Motion  (EOM)  program  used  to 
calculate  the  hydrodynamic  open  water  forces. 
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Q    ********************************************************************** 

SUBROUTINE  TRJ 

q  *************************************************^/g^***************** 

IMPLICIT  NONE 

REAL*8CTRN(6,6),RHS(6),EDOT(6),E(6),PDOT(6) 

REAL*8  P(6),DT,TA,PDOTO,PDDOTO 

REAL*8DS4,DC4,DS5,DC5,DS6,DC6,DT5,DS5DC6,DS5DS6 

REAL*8  PDDOT(6),RHSTOT(6),RHSC(6) 

REAL*8A(6,6),WORK(6,6),ACOF(6),BCOF(6) 

REAL*8  ACP(3),BCP(3),RHSFM(6),RTD,DTR 

REAL*8  XC,YC,ZC,XMTN,XMAX,YMIN,YMAX,ZMrN,ZMAX,D2,HRAD 

REAL*8  VRAD,UHOST,DIA,RHOC,ATLMAX,VLEN 

REAL* 8  FX,FY,FZ,MX,MY,MZ,AED0T1,RHS1,PTF 

REAL*8  ATL,VEX(2),TEMP1,TEMP2 

REAL*8  P1A,P1T,P2A,P2T,P3A,P3T 

REAL*8  P4A,P4T,P5A,P5T,P6A,P6T 

REAL*8  PlDOT,P2DOT,P3DOT 

REAL*8  P4DOT,P5DOT,P6DOT 

REAL*8  PlDDOT,P2DDOT,P3DDOT 

REAL*8  P4DDOT,P5DDOT,P6DDOT 

INTEGER*2I,J,JJ,K,JJJ,n,KODE,JR,JJR,JP,JJP,JYA,JJYA 

INTEGER*2  JX,JJX,JY,JJY,JZ,JJZ 

CHARACTER*  15  LABEL 

REAL*  8  DD(  1800) 

COMMON  /COMDD/  DD 

COMMON/FGEOM/  XC,YC,ZC,XMIN,XMAX,YMIN,YMAX,ZMIN,ZMAX, 
1  D2,HRAD,VRAD,UHOST,DIA,RHOC,ATLMAX,VLEN 
1  ACP,ACOF,BCP,BCOF 

SAVEJ,JJ,JR,JJR,JP,JJP,JYA,JJYA,JX,JJX,JY,JJY,JZ,JJZ 

DATA  J,JJ,JR,JJR,JP,JJP,JYA,JJYA  /9,9,2,0,2,0,2,0/ 

DATA  JX,JJX,JY,JJY,JZ,JJZ  /2,0,2,0,2,0/ 

EQUIVALENCE  (DD(  1),    DT) 
EQUIVALENCE  (DD(  11),   JJJ),(DD(  4),   RTD) 
EQUIVALENCE  (DD(  10),   TA),(DD(  5),   DTR) 
EQUIVALENCE  (DD(  31),EDOT(l)) 
EQUTVALENCE  (DD(  55),   E(  1 )) 
EQUIVALENCE  (DD(  61),   P(  1)) 
EQUIVALENCE  (DD(  71),  RHS(  1)) 
EQUIVALENCE  (DD(  77),PDOT(l)) 
EQUTVALENCE  (DD(107),A(1,1)) 
EQUIVALENCE  (DD(83),PDDOT(l)) 
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DIMENSION  P1A(10),P1T(10),P2A(10),P2T(10),P3A(10),P3T(10) 
DIMENSION  P4A(  10),P4T(  10),P5A(  10),P5T(  10),P6A(  10),P6T(  10) 

IF(TA.EQ.  0)THEN 

C  READ  IN  THE  GEOMETRY  CONSTRAINTS 

WRITE(6,*) '  OPENING  AND  READING  GEOMETRY  LIMITATION  FILE' 

OPEN(  1  ^FILE^GEOMLEV^STATUS^OLD') 

READ(11,*)XC,YC,ZC 

READ(  1 1,*)  XMTN,YMrN,ZMTN 

READ(1 1,*)  XMAX,YMAX,ZMAX 

READ(1 1,*)  D2,HRAD,VRAD,VLEN 

READ(1 1,*)  UHOST,DIA,RHOC,ATLMAX 

CLOSE  (11) 

C  READ  IN  FORCE  MODEL  COEFFICIENTS 

WRITE(6,*)  '  OPENING  AND  READING  LEAST  SQUARES  COEFFICIENT 
FILE' 

OPEN(  12,FILE=,Z2COEF.EXP\  STATUS-'OLD') 

READ  (12,*)  (ACP(H),  H  =  1,3) 

READ  (12,*)  (BCP(H),  H  =  1,3) 

READ  (12,*)  (ACOF(H),  H  =  1,6) 

READ  (12,*)  (BCOF(H),  H  =  1,6) 

CLOSE  (12) 

C  READ  IN  TRAJECTORY  MODIFICATIONS 

WRITE(6,*) '  OPENING  AND  READING  TRAJECTORY  COMMAND  FILE' 
OPEN^^FILE^TRAJ.CM,  STATUS='OLD') 
READ  (14,*)  (VEX(H),  H  =  1,2) 
READ  (14,*) 
READ  (14,*) 

DO  20  I-  1,10 

READ(  14,*)  P1A(I),P1T(I),P2A(I),P2T(I),P3A(I),P3T(I) 
20    CONTINUE 
READ(14,*) 
READ(14,*) 

DO  30  1=  1,10 

READ(  14,*)  P4A(I),P4T(I),P5 A(I),P5T(I),P6A(I),P6T(I) 

P4A(I)  =  P4A(I)  *  DTR 

P5A(I)  =  P5A(I)  *  DTR 

P6A(I)  =  P6A(I)  *  DTR 
30    CONTINUE 
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CLOSE  (14) 

C  COMPLETE  MODEL  ENTRY 
ELSE 
ENDIF 

C  USE  STANDARD  INITIAL  VELOCITY  PROFILE 

C  MODIFY  TRAJECTORY  BY  SCALING  OF  U,V,  AND  W  RELATIVE 

C  TO  THE  HOST  FOLLOWED  BY  ROLL,PITCFLYAW  RATE  COMMANDS 

IF  (TA  .LT.  PIT(JX))  THEN 

JJX  =  JX  -  1 

P1DDOT  =  (PIA(JX)  -  P1A(JJX))/(P1T(JX)  -  PIT(JJX)) 

PI  DOT  =  PI  DOT  +  PlDDOT*DT 

ELSE 

P1DDOT  =  (PIA(JX)  -  P1A(JJX))/(P1T(JX)  -  PIT(JJX)) 

P1DOT  =  P1DOT  +  PlDDOT*DT 

JX  =  JX+1 

ENDIF 

C       WRITE(6,*)  TA,PlDOT,PlDDOT 

IF  (TA  .LT.  P2T(JY))  THEN 

JJY  =  JY  -  1 

P2DDOT  =  (P2A(JY)  -  P2A(JJY))/(P2T(JY)  -  P2T(JJY)) 

P2DOT  =  P2DOT  +  P2DDOT*DT 

ELSE 

P2DDOT  =  (P2A(JY)  -  P2A(JJY))/(P2T(JY)  -  P2T(JJY)) 

P2DOT  =  P2DOT  +  P2DDOT*DT 

JY  =  JY+  1 

ENDIF 

IF  (TA  LT.  P3T(JZ))  THEN 

JJZ  =  JZ  -  1 

P3DDOT  =  (P3A(JZ)  -  P3A(JJZ))/(P3T(JZ)  -  P3T(JJZ)) 

P3DOT  =  P3DOT  +  P3DDOT*DT 

ELSE 

P3DDOT  =  (P3A(JZ)  -  P3A(JJZ))/(P3T(JZ)  -  P3T(JJZ)) 

P3DOT  =  P3DOT  +  P3DDOT*DT 

JZ  =  JZ  +  1 

ENDIF 

IF  (TA  LT.  P4T(JR))  THEN 

J JR  =  JR  -  1 

P4DDOT  =  (P4A(JR)  -  P4A(JJR))/(P4T(JR)  -  P4T(JJR)) 
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P4DOT  =  P4DOT  +  P4DDOT*DT 

ELSE 

P4DDOT  =  (P4A(JR)  -  P4A(JJR))/(P4T(JR)  -  P4T(JJR)) 

P4D0T  =  P4DOT  +  P4DDOT*DT 

JR  =  JR  +  1 

ENDIF 

IF  (TA  LT.  P5T(JP))  THEN 

JJP  =  JP  -  1 

P5DDOT  =  (P5A(JP)  -  P5A(JJP))/(P5T(JP)  -  P5T(JJP)) 

P5DOT  =  P5DOT  +  P5DDOT*DT 

ELSE 

P5DDOT  =  (P5A(JP)  -  P5A(JJP))/(P5T(JP)  -  P5T(JJP)) 

P5DOT  =  P5DOT  +  P5DDOT*DT 

JP  =  JP  +  1 

ENDIF 

IF  (TA  LT.  P6T(JYA))  THEN 

JJYA  =  JYA  -  1 

P6DD0T  =  (P6A(JYA)  -  P6A(JJYA))/(P6T(JYA)  -  P6T(JJYA)) 

P6DOT  =  P6DOT  +  P6DDOT*DT 

ELSE 

P6DDOT  =  (P6A(JYA)  -  P6A(JJYA))/(P6T(JYA)  -  P6T(JJYA)) 

P6DOT  =  P6DOT  +  P6DD0T*DT 

JYA  =  JYA+  1 

ENDIF 

C      ESTABLISH  THE  BODY  POSITION,  VELOCITY,  AND  ACCELERATION 
C      OVER  TIME  IN  INERTIAL  REFERENCE  FRAME  COORDINATES 
C      RELATIVE  TO  THE  HOST 

PDOTO  =  VEX(1)*(1-  DEXP(-(VEX(2))*(TA**2))) 
PDOT(l)  =  PDOTO*DCOS(P(6))  +  P1DOT 
PDOT(2)  =  PDOTO*DSIN(P(6))  +  P2DOT 
PDOT(3)  =  PDOTO*DTAN(P(5))  +  P3DOT 
PDOT(4)  =  P4DOT 
PDOT(5)  =  P5DOT 
PDOT(6)  =  P6DOT 

PDDOTO  =  VEX(1)*(2.*VEX(2)*TA*DEXP(-(VEX(2))*(TA**2))) 

PDDOT(l)  =  PDDOTO*DCOS(P(6))  +  P1DDOT 

PDDOT(2)  =  PDDOTO*DSIN(P(6))  +  P2DDOT 

PDDOT(3)  =  PDDOTO*DTAN(P(5))  +  P3DDOT 

PDDOT(4)  =  P4DDOT 

PDDOT(5)  =  P5DDOT 
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PDDOT(6)  =  P6DDOT 

P(l)  =  P(l)  +  PDOT(l)*DT  +  PDDOT(l)*DT**2 
P(2)  =  P(2)  +  PDOT(2)*DT  +  PDDOT(2)*DT**2 
P(3)  =  P(3)  +  PDOT(3)*DT  +  PDDOT(3)*DT**2 
P(4)  =  P(4)  +  PDOT(4)*DT  +  PDDOT(4)*DT**2 
P(5)  =  P(5)  +  PDOT(5)*DT  +  PDDOT(5)*DT**2 
P(6)  -  P(6)  +  PDOT(6)*DT  +  PDDOT(6)*DT**2 

C  CALCULATE  COORDINATE  TRANSFORM  FROM  INERTIA  FRAME  (HOST) 
C  TO  VEHICLE  FRAME 

C  TRANSPOSE  OF  THE  SNAME  1950  CONVENTION  (LAW  ABKOWICZ) 
C  ORDER  IS  YAW,  PITCH,  ROLL 

DS4  =  DSIN(P(4)) 
DC4  =  DCOS(P(4)) 
DS5  -  DSIN(P(5)) 
DC5  =  DCOS(P(5)) 
DS6  =  DSIN(P(6)) 
DC6  =  DCOS(P(6)) 


CTRN(1,L 

>=  DC5*DC6 

CTRN(1,2; 

>=  DC5*DS6 

CTRN(1,3; 

>  =  -DS5 

CTRN(2,r 

>  =  -DS6*DC4  +  DS4*DS5*DC6 

CTRN(2,2; 

)  =  DC4*DC6  +  DS4*DS5*DS6 

CTRN(2,3; 

)=  DS4*DC5 

CTRN(3,L 

»=  DS4*DS6 +DC4*DC6*DS5 

CTRN(3,2; 

)  =  -DS4*DC6  +  DC4*DS5*DS6 

CTRN(3,3; 

)=  DC4*DC5 

CTRN(4,4; 

1=  1.0 

CTRN(4,5; 

1=  0.0 

CTRN(4,6; 

1  =  -DS5 

CTRN(5,4; 

>=  0.0 

CTRN(5,5; 

>=  DC4 

CTRN(5,6; 

}=  DS4*DC5 

CTRN(6,4; 

)=  0.0 

CTRN(6,5; 

)  =  -DS4 

CTRN(6,6 

)=  DC5*DC4 

C  XFORM  INERTIA  RELATIVE  VELOCITIES  AND  ACCELERATIONS 

C  TO  VEHICLE  RELATTVE 

C      CTRN(6X6)*PDOT(6Xl)  =  E(6Xl) 

C      CTRN(6X6)*PDD0T(6X1)  =  ED0T(6X1) 


68 


C  INITIALIZE  THE  E  AND  EDOT  MATRICE  TO  ZERO  AT  EACH  TIME  STEP 

DO  110  1=1,6 
E(I)  =  0.0 
EDOT(I)  =  0.0 
110  CONTINUE 

DO  120  1=1,6 

DO  130  K=l,6 

E(I)  =  E(I)  +  PDOT(K)*CTRN(I,K) 

EDOT(I)  =  EDOT(I)  +  PDDOT(K)*CTRN(I,K) 
130  CONTINUE 
120  CONTINUE 

C  CONVERT  BODY  VELOCITIES  TO  ADD  IN  HOST  VESSEL  SPEED  EFFECTS 
C  HOST  SPEED  IN  (-)  P(l)  DIRECTION  ONLY.  ADD  AS  A  POSITIVE  SCALAR 

E(l)  =  E(l)  +  UHOST*DCOS(P(6))*DCOS(P(5)) 

E(2)  =  E(2)  -  UHOST*DSIN(P(6)) 

E(3)  =  E(3)  +  UHOST*DCOS(P(6))*DSIN(P(5)) 

C  SOLVE  FOR  THE  RIGHT  HAND  SIDE  VECTOR  WITH  TOTAL  FORCE 
C  EDOT  (6X1)  *  LHS  (6X6)  =  RHS  (6X1) 

C  INITIALIZE  THE  RHSTOT  TO  ZERO  AT  EACH  TIME  STEP 
DO  140  1=1,6 
RHSTOT(I)  =  0.0 
140  CONTINUE 

IF  (JJ  EQ.  9)  THEN 

WRITE(25,178)  (acp(I),I=l,3) 

LABEL  =  'CTRK 

CALL  MPRINT1(CTRN,6,6,9,LABEL) 

LABEL  =  'LHS' 

CALL  MPRINT1(A,6,6,9,LABEL) 

ENDIF 

JJ  =  JJ  +  1 

DO  200  1=1,6 

DO  220  K=  1,6 

RHSTOT(I)  =  RHSTOT(I)  +A(I,K)*EDOT(K) 
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220  CONTINUE 
200  CONTINUE 

AEDOTl=RHSTOT(l) 

C  CALCULATE  THE  RIGHT  HAND  SIDE  OF  THE  HYDRODYNAMIC 

EQUATIONS 

C  WITH  NO  CONTROL  FORCE  OR  INTERACTION  FORCE 

C  SET  CROSSFLOW  SWITCH  JJJ  TO  1 

JJJ=  1 

CALL  HYDRO 

C  CALCULATE  THE  INTERACTION  FORCES  IN  FORCE  MODEL  ROUTINE 

CALL  FM  (FX,FY,FZ,MX,MY,MZ,AEDOTl,PTF,ATL,KODE) 

C  CALCULATE  THE  CONTROL  FORCES  TO  EFFECT  THE 
C  CHOSEN  TRAJECTORY.  ADD  IN  HYDRO  FORCE 

RHSFM(1)  =  FX 
RHSFM(2)  =  FY 
RHSFM(3)  =  FZ 
RHSFM(4)  =  MX 
RHSFM(5)  =  MY 
RHSFM(6)  =  MZ 

DO  240  1=1,6 

RHSC(I)  =  RHSTOT(I)  -  RHSFM(I) 
240  CONTINUE 

C  PRINT  OUTPUT  EVERY  10  DELTA  T 

J  =  J+1 

JJF(J.EQ.  10.)  THEN 

J  =  0 

WRITE(2 1,176)  TA,  (RHSTOT(I),I=l,6) 
WRITE(22,178)  TA,  KODE,  (RHS(I),I=1,6) 
WRITE(23,178)  TA,  KODE,  (RHSFM(I),I=1,6) 
WRITE(24,177)  TA,  KODE,  (RHSC(I),I=1,6) 
WWTE(25,176)TA,(P(I),I=1,3),(RTD*P(I),I=4,6) 
WRITE(26, 176)  TA,(PDOT(I),I=l,3),(RTD*PDOT(I),I=4,6) 
WRITE(27,176)TA,(PDDOT(I),I=l,3),(RTD*PDDOT(I),I=4,6) 
WRITE(28,176)TA,(E(I),I-l,3),(RTD*E(I),I-4,6) 
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WRITE(29, 176)  TA,(EDOT(I),I=l,3),(RTD*EDOT(I),I=4,6) 
176  F0RMAT(F5.1,1X,6(E11.4E3,1X)) 
177FORMAT(F6.1,lX,I2,lX,6(E11.3E3,lX)) 
178F0RMAT(F6.1,1X,I2,1X,6(E9.3,1X)) 

ENDIF 


RETURN 
END 
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/-<  ********************************************************************** 

SUBROUTINE  FM  (FX,FY,FZ,MX,MY,MZ, 
1  AEDOTl,PTF,ATL,KODE) 

C  SUBTOUTINE  FORCE  MODEL  --  VEHICLE  SPECIFIC  FORCE  AND  MOMENT 

C  MODELS  ARE  PLACED  HERE 

C 

C  CENTER  OF  X,Y,Z,  COORDINATE  SYSTEM  IS  ON  HOST  AXIS 

C  IN  THIS  VERSION,  IT  IS  PRESUMED  THAT  A  NON-CRASH  TRAJECTORY 

C  SPECIFIED 

C  MODEL  IS  VALID  FOR  THE  STARBOARD  HALF  PLANE 

C  (ALPHA  CALCULATION) 

REAL*8  MX,MY,MZ,MYC,MZC,MS,MTOT 

REAL*8  P(6),E(6),XC,YC,ZC,XMrN,YMTN,ZMIN 

REAL*8  XMAX,YMAX,ZMAX,D2,HRAD,VRAD,UHOST,DIA 

REAL*8RHOC,ATLMAX,VLEN,ACP(3),ACOF(6),BCP(3) 

REAL*8BCOF(6),FX,FY,FZ,DOGS,PTF,PTFS,PTFC 

REAL*8UZ1LF,UZ1LM,UZ1VF,UZ1VM,FACC,FACS,GAPS 

REAL*8FYC,FZC,GAP,DOG,PT,SG,FTOT,POSRAD 

REAL*8UZ2LF,UZ2LM,X,Y,Z,ATL,ATLINV,AXDRAG 

REAL*8  U,  V,W,PHI,THETA,PSI,ALPHA,FAREA,C  1,C2,C4,C5 

REAL*8  Dl,RAT,XS,YS,ZS,AEDOTl,RHSl,TA 

REAL*8  RHS(6) 

REAL*8DD(1800) 

REAL*4  FOO,F001,F002,PTF2 

INTEGER  KODE 
COMMON  /COMDD/  DD 

COMMON  /FGEOM/  XC,YC,ZC,XMTN,XMAX,YMIN,YMAX,ZMrN,ZMAX, 
1  D2,HRAD, 
WAD,UHOST,DIA,RHOC,ATLMAX,VLEN,ACP,ACOF,BCP,BCOF 
EQUIVALENCE  (DD(  55),    E(1)),(DD(  61),    P(l)) 
EQUIVALENCE  (DD(  10),      TA),(DD(  71),  RHS(l)) 


PI  =  3. 14159 
DPR=  180.0/pi 
X  =  P(1) 

Y  =  P(2) 
Z  =  P(3) 
U  =  E(1) 

V  =  E(2) 
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W  =  E(3) 

PHI  =  P(4) 

THETA  =  P(5) 

PSI  =  P(6) 

ALPHA  =  ATAN2(Z,Y) 

FX  =  0.0 

FY -0.0 

FZ  =  0.0 

MX  =  0.0 

MY  =  0.0 

MZ  =  0.0 

FAREA  =  PI/4.*DIA**2 

CI  =  5*RHOC*(UHOST**2)*FAREA 

C2  =  C1*VLEN 

C4  =  5*(RHOC*U**2)*FAREA 

C5  =  C4*VLEN 

C6  =  C5/12. 

C  IN  THIS  SUBROUTINE,  THE  HOST  IS  TREATED  AS  A  CYLINDER  OF  RADIUS 

HRAD. 

C  POSRAD  IS  RADIUS  FROM  HOST  AXIS  TO  SMALL  VEHICLE  ORIGIN. 

POSRAD  =  SQRT(Y*Y  +  Z*Z) 

Dl  =  SQRT(  (X-XC)**2  +  (Y-YC)**2  +  (Z-ZC)**2  ) 

IF  ((Z  LE.  ZC)  AND.  (Y  LE.  YMAX)  AND.  (Y  GE.  YMIN)  AND. 
1      (X  LE.  XMAX)  AND.  (X  GE.  XMIN))  THEN 

C  VEHICLE  IS  IN  THE  LAUNCHWAY  -  ZONE  1  * 

ATL  =  ATLMAX-(SQRT(  (X-XC)**2  +  (Y-YC)**2  +  (Z-ZC)**2))/VLEN 

ATLINV  =  ATL  -  ATLMAX 

AXDRAG  =  (.637341  +.853415*(ATLINV)**3)*C1 

IF  (AXDRAG  It.  0.0)  AXDRAG  =  0.0 

FTOT  =  AEDOT1  +  AXDRAG 

PTF  =  .5*(FTOT)/Cl 

FX  =  -AXDRAG 

UZ1LF  =  ((-2.04*ATL**3+5.59*ATL**2-4.49*ATL+1.13)*PTF  + 
1(1.40657-.55068*ATL)) 

FY  =  UZ1LF*C1 
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UZ1LM  =  -((.  153*ATL**3  -.198*ATL**2  -  .516*ATL  +.508)*PTF) 
MZ  =  UZ1LM*C2 
KODE  =  1 


ELSEBF  (  (Z  LE.  ZC)  AND.  (Y  GE.  YMAX)  AND.  (X  GE.  XMAX) 
1         AND.  (Dl  XT.  D2))  THEN 

C  VEHICLE  IS  IN  TRANSITION  ZONE  -  HEMI-SPHERE  ZONE  2      * 

C  * 

IF  (Dl  XT.  0.001)  Dl=  0.001 

RAT  =  D2/D1 

XS  =  XC  +  RAT*(X-XC) 

YS  =  YC  +  RAT*(Y-YC) 

ZS  =  ZC  +  RAT*(Z-ZC) 

FACC  =  (D2-D1)/D2 

FACS  =  D1/D2 

GAPS  -  SQRT(YS*YS  +  ZS*ZS)  -  HRAD  -  VRAD 

IF  (GAPS  .It.  0.001)  GAPS  -  0.001 

DOGS  =  DIA/GAPS 

ATL=  ATLMAX 

ATLINV  =  0.0 

AXDRAG  =  (.637341)*C1 

FXC  =  AEDOT1  +  AXDRAG 

C  AT  LAUNCHER  EXIT 

PTFC  =  5*(FXC)/C1 


FYC  =  UZ1LF*C1 

UZ1LF  =  ((-2.04*ATL**3+5.59*ATL**2-4.49*ATL+1.13)*PTFC  + 
1     (1.40657-.55068*ATL)) 

UZ1LM  =  -((.153*ATL**3  -.198*ATL**2  -  .516*ATL +.508)*PTFC) 

MZC  =  UZ1LM*C2 
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C  PTF  EFFECT  AT  SPHERE  INTERSECTION  POINT 

FXS  =  -RHS(l)  +  AEDOT1 

PTFS  -  5*(FXS)/C4 

UZ2LF  =  -(0.004*DOGS*PTFS) 

UZ2LM  =  (4. 1 3  743  *DEXP(- 1 56. 1 97*(  1  ./DOGS)))*PTFS 

Z2LMPTF  =  UZ2LM*C5 

C  ADD  WALL  RELATIVE  PITCH  AND  GAP  EFFECT 

C  APPLY  EFFECTIVE  PITCH  ANGLE  DUE  TO  LAUNCHWAY  GEOMETRY 

PTEFF  =  ((-THETA*SIN(ALPHA)  +  PSI*COS(ALPHA))  *  DPR) 
PT  =  PTEFF*  (XS-XC)/(2.5*VLEN) 
CDBG      WRITE(6,*)  PTEFF,PT 

SG=SQRT(DOGS) 

UZ2LF  =  UZ2LF 
1     -((ACP(l)*DEXP(ACP(2)*(l./DOGS)**ACP(3)))+ 
1     (ACOF(l)*PT*SG+ACOF(2)*PT*SG*SG+ACOF(3)*PT*SG**3  + 
1     ACOF(4)*PT*PT*SG+ACOF(5)*PT*PT*SG*SG+ACOF(6)*PT*PT*SG**3)) 

FS  =  UZ2LF*C4 

UZ2LM  = 
1     (BCOF(l)*PT*SG+BCOF(2)*PT*SG*SG+BCOF(3)*PT*SG**3  + 
1     BCOF(4)*PT*PT*SG+BCOF(5)*PT*PT*SG*SG+BCOF(6)*PT*PT*SG**3) 

MS  =  UZ2LM*C5+Z2LMPTF 

FX  =  FACC*(-FXC)  +  FACS*(-FXS) 

FY  =  FACC*FYC  -  FACS*FS*COS(ALPHA)  +  FACS*RHS(2) 
FZ  =  FACC*FZC  -  FACS*FS*SIN(ALPHA)  +  FACS*RHS(3) 

MY  =  FACC*MYC  -  FACS*MS*SIN(ALPHA)  +  FACS*RHS(5) 

MZ  =  FACC*MZC  +  FACS*MS*COS(ALPHA)  +  FACS*RHS(6) 

PTF  -  FACC*PTFC  +  FACS*PTFS 

KODE  =  3 
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c 

ELSE 

C    ZONE  IS  ZONE  2  * 

GAP  =  SQRT(Y*Y  +  Z*Z)  -  HRAD  -  VRAD 

IF  (GAP  .LT.  0.001)  GAP  =  0.001 

DOG  =  DIA/GAP 

FX  =  -RHS(l)  +  AEDOTl 

C  PTF  EFFECT 

PTF  =  5*(-RHS(l)  +  AEDOTl)/C4 

UZ2LF  =  -(0.004*DOG*PTF) 

UZ2LM  =  (4. 13743*DEXP(-156. 197*(l./DOG)))*PTF 

Z2LMPTF  =  UZ2LM*C5 

C  ADD  WALL  RELATIVE  PITCH  AND  GAP  EFFECT 

PT  =  (-THETA*SIN( ALPHA)  +  PSI*COS(ALPHA))  *  DPR 
SG-SQRT(DOG) 

UZ2LF  =  UZ2LF 
1     -((ACP(  1)*DEXP(ACP(2)*(  l./DOG)**ACP(3)))+ 
1     (ACOF(l)*PT*SG+ACOF(2)*PT*SG*SG+ACOF(3)*PT*SG**3  + 
1     ACOF(4)*PT*PT*  SG+ACOF(5  )*PT*PT*  SG*  SG+ACOF(6)*PT*PT*  SG*  *3  )) 

FTOT  -  UZ2LF*C4 

UZ2LM  = 
1     (BCOF(l)*PT*SG+BCOF(2)*PT*SG*SG+BCOF(3)*PT*SG**3  + 
1     BCOF(4)*PT*PT*SG+BCOF(5)*PT*PT*SG*SG+BCOF(6)*PT*PT*SG**3) 

MTOT  =  UZ2LM*C5  +  Z2LMPTF 

FX=  -FX 

FY  =  -  FTOT*COS(ALPHA)  +  RHS(2) 

FZ  =  -  FTOT*SIN(ALPHA)  +  RHS(3) 


76 


MY  =  -  MTOT*SIN(ALPHA)  +  RHS(5) 

MZ  =     MTOT*COS(ALPHA)  +  RHS(6) 

KODE  =  2 

ENDIF 

RETURN 

END 
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SUBROUTINE  HYDRO 

/-<  ********************************************************************** 

C  DTRC  2510  HYDRODYNAMIC  COEFFICIENTS 

C  NOW  INCLUDES  ZERO  FORWARD  SPEED  COEFFICIENTS 

C 

IMPLICIT  NONE 

INTEGER*2  JJIILJ 

INTEGER*2  IDTAR 

REAL* 8  TEMP  1,TEMP2,TEMP3, TEMP 

REAL*8UPFRAC,LOFRAC,URACER,URACER0,URACES,URACES0 

REAL*8DRE3,DSPE3,DSSE3,DRE,DSPE,DSSE 

REAL*8  DRHS  1,DRHS2,DRHS3,DRHS4,DRHS5,DRHS6 

REAL*8  ZRHS1,ZRHS2,ZRHS3,ZRHS4,ZRHS5,ZRHS6 

REAL*8  BRHS  1,BRHS2,BRHS3,BRHS4,BRHS5,BRHS6 
C 

REAL*  8  DD(  1800) 

COMMON  /COMDD/  DD 
C 

REAL*8URUR,USUS,DELRCR,DELRCS,ACOF(8),BCOF(8),acp(4),bcp(4) 

REAL*8DSP,DSS,DR,DB,UC,DT,RHS(6),E(6),P(6),QPROP,TA,YG,YB 

REAL*8PHISTERN,VATS,WATS,ETA,CETA,ETAM1,BETAS 

REAL*8  KBETA,WT,BOUY 

REAL*8  ARAPPHI2,ARAPPffl4,CROSSLIM,SINTHETA 

REAL*8  COSTHETA,SINPHLCOSPHI 

REAL*8  UU,W,WW,PP,QQ,RR,PR,PQ,QR,UV,UW,UP 

REAL*  8  UQ,UR,VATTAIL,UUDB 

REAL*8  VP,VR,WP,WQ,SQRWWW,VSQRWWW 

REAL*8  WSQRWWW,CROSS(6),VORT(6) 

REAL* 8  DRDR,DSDS,UUDS,UUDR,UABW,UABV,WABW 

REAL* 8  QABQ,VABV,RABR,FXP 

REAL*  8  PABP,DPDP,UUDP,UT,BETA, ALPHA, WR 

REAL*  8  ARRAY(6),WTINH20,DBDB 
C 

REAL*8  XUDOT,XQQ,XRP,XVR,XWQ,XW,XDRDR,XDSSDSS 

REAL* 8  XDSPDSP,XDBDB 

REAL*8XW,XQ,XWDSP,XWDSS,XQDSS,XQDSP,XVDR,XRDR,XWW,XRR 

REAL*8  XDRE2,XDSSE2,XDSPE2 

REAL*  8  YVDOT,YRDOT,YRU,YVU,YVABVU,YW 

REAL*  8  YDR,YDRN,YRABR,YVABV 

REAL*  8  YPU, YWP, YPABP, YPQ, YSTAR, YVAB VR, YHVAB V, YDRE, YDRE3 

REAL*  8  ZWDOT,ZQDOT,ZQU,ZWU,ZABWU,ZWW 

REAL*8  ZDSS,ZDSSN,ZQABQ,ZWABW 

REAL*8  ZVP,ZSTAR,ZPR,ZDSP,ZDSPN,ZHWABW,ZDB 

REAL*8  ZDSSE,ZDSPE,ZDSPE3,ZDSSE3 
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REAL* 8  KPDOT,KVDOT,KRDOT,KQR,KPU,KRU,KWP,KVU,K4S 

REAL*8  K8S,KI,KDSP,KDSS 

REAL*8  KPABP,KWR,KDR,KDRN,KSTAR 

REAL*8KDRE,KDRE3,KDSSE,KDSPE,KDSSE3,KDSPE3 

REAL*8  MQDOT,MWDOT,MQU,MWU,MWABWR,MABWU,MWW,MDSS 

REAL*8  MDSSN,MQABQ,MWABW 

REAL*8  MRP,MSTAR,MVP,MDSP,MDSPN,MHWABW,MDB 

REAL*8  MDSSE,MDSPE,MDSSE3,MDSPE3 

REAL*  8  NRDOT,NVDOT,NRU,NVU,NVABVR,NABVU,NW,NDR 

REAL*  8  NDRN,NVABV,NRABR 

REAL* 8  NWP,NPQ,NSTAR,NPU,NHVABV,NDRE,NDRE3 

REAL*8IYY,IZZ,LXX,SLEN,MASS,ZB,ZG,XB,XG,XSTERN 

REAL*8XFORWARD,KTAINF,UAINF,IXZ,LXY,IYZ,YPDOT,NPDOT 

EQUIVALENCE  (DD(   1),      DT),  (DD(  10),      TA) 

EQUIVALENCE  (DD(  11),     JJJ),  (DD(  16),    BETA) 

EQUIVALENCE  (DD(  15),  ALPHA) 

EQUIVALENCE  (DD(  17),     FXP  ),(DD(  18),    QPROP) 

EQUIVALENCE  (DD(  19),     ETA  ),(DD(  20),     CETA) 

EQUIVALENCE  (DD(  55),     E(1)),(DD(  61),     P(l)) 

EQUIVALENCE  (DD(  67),      DSS),(DD(  68),      DSP) 

EQUIVALENCE  (DD(  69),       DR),(DD(  71),  RHS(l)) 

EQUIVALENCE  (DD(  1 79),       DB) 

EQUIVALENCE  (DD(260),CROSS(1)),(DD(266),  VORT(l)) 

EQUIVALENCE  (DD(272),VATTADL) 

EQUIVALENCE  (DD(301),XUDOT    ),(DD(303),XW       ) 
EQUIVALENCE  (DD(305),XQ       ),(DD(307),XQQ      ) 
EQUIVALENCE  (DD(309),XRR      ),(DD(311),XRP      ) 
EQUIVALENCE  (DD(3 1 3 ),XVR      ),(DD(3 1 5 ),XWQ      ) 
EQUIVALENCE  (DD(317),XW      ),(DD(319),XWW      ) 
EQUIVALENCE  (DD(321),XDBDB    ),(DD(323),XDRDR    ) 
EQUIVALENCE  (DD(325),XDSSDSS  ),(DD(327),XDSPDSP  ) 
EQUIVALENCE  (DD(329),XDRE2    ),(DD(331),XDSSE2  ) 
EQUIVALENCE  (DD(333),XDSPE2  ) 

EQUIVALENCE  (DD(1 703 ),XWDSP  ),(DD(1705),XWDSS  ) 
EQUIVALENCE  (DD(1707),XQDSP  ),(DD(1709),XQDSS  ) 
EQUTVALENCE(DD(1711),XVDR    ),(DD(1713),XRDR    ) 

EQUIVALENCE  (DD(3 38), YVDOT    ),(DD(340),YPDOT    ) 
EQUIVALENCE  (DD(342),YRDOT    ),(DD(344),YPABP    ) 
EQUIVALENCE  (DD(346),YPQ      ),(DD(348),YRU      ) 
EQUIVALENCE  (DD(350),YPU      ),(DD(352),YVU      ) 
EQUIVALENCE  (DD(3  54), YVABVR  ),(DD(356),YVABVU  ) 
EQUIVALENCE  (DD(3 5 8), YW      ),(DD(360),YDR      ) 
EQUIVALENCE  (DD(362),YDRN     ),(DD(364),YRABR    ) 
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EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 


DD(366),YVABV    ),(DD(368),YWP      ) 
DD(370),YSTAR    ),(DD(372),YHVABV  ) 
DD(374),YDRE3    ),(DD(376),YDRE     ) 

DD(390),ZWDOT    ),(DD(392),ZQDOT    ) 
DD(394),ZQU      ),(DD(396),ZVP      ) 
DD(398),ZSTAR    ),(DD(400),ZWU      ) 
DD(402),ZABWU    ),(DD(404),ZWW      ) 
DD(406),ZDB      ),(DD(408),ZDSS     ) 
DD(410),ZDSSN    ),(DD(412),ZDSP     ) 
DD(414),ZDSPN    ),(DD(416),ZQABQ    ) 
DD(418),ZWABW    ),(DD(420),ZHWABW  ) 
DD(  422),ZDSSE3  ),(DD(  424),ZDSPE3  ) 
DD(  426),ZDSSE  ),(DD(  428),ZDSPE  ) 
DD(1717),ZPR      ) 

DD(440),KVDOT    ),(DD(442),KPDOT    ) 
DD(444),KRDOT    ),(DD(446),KQR      ) 
DD(448),KPABP    ),(DD(450),KPU      ) 
DD(452),KRU      ),(DD(454),KWP      ) 
DD(456),KSTAR    ),(DD(458),KVU      ) 
DD(460),KI       ),(DD(462),K4S      ) 
DD(464),K8S      ),(DD(466),KDR      ) 
DD(468),KDRN     ),(DD(470),KDSS     ) 
DD(472),KDSP     ),(DD(474),KDRE     ) 
DD(476),KDRE3    ),(DD(478),KDSSE    ) 
DD(480),KDSPE    ),(E>D(482),KDSSE3  ) 
DD(484),KDSPE3  ) 
DD(1721),KWR     ) 

DD(494),MWDOT    ),(DD(496),MQDOT    ) 
DD(498),MRP      ),(DD(500),MQU      ) 
DD(502),MWU      ),(DD(504),MSTAR    ) 
DD(506),MWABWR  ),(DD(508),MABWU    ) 
DD(510),MWW      ),(DD(512),MDSS     ) 
DD(514),MDSSN    ),(DD(516),MDSP     ) 
DD(518),MDSPN    ),(DD(520),MDB      ) 
DD(522),MQABQ    ),(DD(524),MWABW    ) 
DD(  526),MHWABW  ),(DD(  528),MDSSE3  ) 
DD(  530),MDSPE3  ),(DD(  532),MDSSE  ) 
DD(  534),MDSPE  ) 
DD(1723),MVP      ) 

DD(546),NVDOT    ),(DD(548),NPDOT    ) 
DD(550),NRDOT    ),(DD(552),NPQ      ) 
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EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 

EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 
EQUIVALENCE 


(DD(554) 
(DD(558) 
(DD(562) 
(DD(566) 
(DD(570) 
(DD(574) 
(DD(578) 
(DD(582) 

(DD(593) 
(DD(598) 
(DD(602) 
(DD(606) 
(DD(610) 
(DD(614) 
(DD(618) 
(DD(622) 
(DD(626) 
(DD(629) 
(DD(630) 
(DD(631) 


NPU      ),(DD(556),NRU      ) 
NSTAR    ),(DD(560),NVU      ) 
NVABVR  ),(DD(564),NABVU 
NW      ),(DD(568),NDR      ) 
NDRN     ),(DD(572),NVABV    ) 
NRABR    ),(DD(576),NWP      ) 
NHVABV  ),(DD(580),NDRE3 
NDRE     ) 


) 


SLEN 

IZZ 

LXZ 

LXY 

WT 

ZB 

XB 


) 


) 


),(DD(596),IYY 
),(DD(600),LXX 
),(DD(604),IYZ      ) 
),(DD(608),MASS     ) 
),(DD(612),BOUY     ) 
),(DD(616),ZG       ) 
),(DD(620),XG       ) 

XSTERN  ),(DD(624),XFORWARD) 

YG       ),(DD(628),YB       ) 

WTINH20) 

CROSSLIM) 

LOFRAC  ),(DD(632),UPFRAC  ) 


(DD(848),KTAINF  ),(DD(849),UAINF    ) 
(DD(850),URACERO  ),(DD(851),URACES0  ) 
(DD(858),URACER  ),(DD(859),URACES  ) 
(DD(866),DRE      ),(DD(867),DSSE     ) 
(DD(868),DSPE     ) 

(DD(  1 75  8),  ARRAY(  1  )),(DD(  1 764),  ACOF(  1 )) 
(DD(1800),botdep),  (DD(1772),BCOF(l)) 
(DD(1780),acp(l)),  (DD(1788),  bcp(l)) 


C  FOR  EFFICIENCY  THE  FOLLOWING  VARIABLES  ARE  CALCULATED 
C 

UU  =  E(1)*E(1) 

W  =  E(2)*E(2) 

WW  =  E(3)*E(3) 

UT  =  DSQRT(UU+W+WW) 

PP  =  E(4)*E(4) 

QQ  =  E(5)*E(5) 

RR  =  E(6)*E(6) 

PR  -  E(4)*E(6) 

PQ  -  E(4)*E(5) 

QR  =  E(5)*E(6) 

UV-E(1)*E(2) 
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UW  =  E(1)*E(3) 

UP  =  E(1)*E(4) 

UQ  =  E(1)*E(5) 

UR  =  E(1)*E(6) 

VP  =  E(2)*E(4) 

VR  =  E(2)*E(6) 

WP  =  E(3)*E(4) 

WQ  =  E(3)*E(5) 

WR  =  E(3)*E(6) 

URUR  =  URACER*URACER 

USUS  =  URACES*URACES 

DELRCR  =  URACER  -  URACERO 

DELRCS  =  URACES  -  URACESO 

VATS  =  E(2)  +  XSTERN*E(6) 

WATS  =  E(3)  -  XSTERN*E(5) 

IF(  DABS(URACER)  .GT.  0. 1)  THEN 
DRE  =DR-DATAN2(  VATS,URACER ) 

ELSE 
DRE  =  DR 

ENDIF 

IF(  DABS(URACER)  .GT.  0.1)  THEN 
DSSE  -  DSS  -  DATAN2(  WATS,URACES  ) 
DSPE  =  DSP  -  DATAN2(  WATS,URACES  ) 
ELSE 

DSSE  =  DSS 
DSPE  -  DSP 

ENDIF 

DRE3    =DRE*DRE*DRE 

DSPE3  =  DSPE*DSPE*DSPE 

DSSE3  =  DSSE*DSSE*DSSE 

DRDR  =  DR*DR 

DBDB  =  DB*DB 

DSDS  =  DSS*DSS 

DPDP  =  DSP*DSP 

UUDR  =  UU*DR 

UUDB  =  UU*DB 

UUDP  =  UU*DSP 

UUDS  =  UU*DSS 

UABW  =  E(1)*ABS(E(3)) 

UABV  =  E(1)*ABS(E(2)) 

WABW  =  E(3)*ABS(E(3)) 

QABQ  =  E(5)*ABS(E(5)) 


82 


VABV  =  E(2)*ABS(E(2)) 
PABP  =  E(4)*ABS(E(4)) 
RABR  =  E(6)*ABS(E(6)) 
SQRWWW    =  DSQRT(W+WW) 
VSQRWWW  =  E(2)  *  SQRWWW 
WSQRWWW  =  E(3)  *  SQRWWW 
WTINH20    =WT-BOUY 
SINPHI     =  DSIN(P(4)) 
COSPHI     =  DC0S(P(4)) 
SINTHETA  =  DSIN(P(5)) 
COSTHETA  =DCOS(P(5)) 
C 

IF(  (DABS(VATS)  .GT.  0.001)  .OR.  (DABS(WATS)  .GT.  0.001))  THEN 

TEMPI  =  VATS*  VATS 

TEMP2  =  WATS*WATS 

TEMP3  =  TEMP1+TEMP2 

ARAPPHI2  =  TEMP3*2.*VATS*WATS/TEMP3 

ARAPPffl4  =  TEMP3*4.*VATS*WATS*(TEMP2  -  TEMP1)/(TEMP3*TEMP3) 

BETAS     =  DATAN2(DSQRT(TEMP3),E(1)) 

PfflSTERN  =  DATAN2(WATS,-VATS) 
ELSE 

ARAPPHI2  =  0.0 

ARAPPHI4  -  0.0 

BETAS     -  0.0 

PfflSTERN  =  0.0 
ENDIF 
C 

IF((DABS(E(3))  .GT.  0.1D-10)  .AND.  (DABS(E(1))  .GT.  0.01))  THEN 

ALPHA  =  DATAN2(E(3),E(1)) 
ELSE 

ALPHA  =  0.0 
ENDIF 
C 

IF(  (DABS(E(2))  .GT.  0.1D-20)  .AND.  (DABS(UT)  .GT.  .01))  THEN 

BETA  =  DASIN(-E(2)/UT) 
ELSE 

BETA  =  0.0 
ENDIF 
C 

C  FOR  EFFICIENCY  CALL  THESE  SUBROUTINES  ONLY  ONCE 
C 

IF(JJJ.EQ.  1)  THEN 

CALL  VORTEX 

CALL  CROSFLOW 

CALL  PROP 
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c 

D0  345IDTAR  =  1,10 
CALL  TOWDARRY(IDTAR) 
345    CONTINUE 

ENDIF 
C 

CETA 

ETAM1  =  ETA  -  1.0 
C 

C  CALCULATE  RIGHT  HAND  SIDES  OF  EQUATIONS  OF  MOTION 
C 

C  AXIAL  FORCES  ZERO  SPEED  MODEL 
C 

ZRHS1  =  0.0 
C 

C  LATERAL  FORCE  ZERO  SPEED  MODEL 
C 

ZRHS2  =  YHVABV*VATS*DABS(VATS) 
C 

C  NORMAL  FORCE  ZERO  SPEED  MODEL 
C 

ZRHS3  =  ZHWABW*WATS*DABS(WATS) 
C 

C  ROLL  MOMENT  ZERO  SPEED  MODEL 
C 

ZRHS4  =  KPABP*PABP 
C 

C  PITCH  MOMENT  ZERO  SPEED  MODEL 
C 

ZRHS5  =  MHWABW*WATS*DABS(WATS) 
C 

C  YAW  MOMENT  ZERO  SPEED  MODEL 
C 

ZRHS6  =  NHVABV*VATS*DABS(VATS) 
C 

C  CALCULATE  DTRC  HYDRODYNAMIC  FORCES 
C  HYDRODYNAMIC  AXIAL  FORCES 
C 

DRHS 1  =  XQQ*QQ  +  XRR*RR  +  XRP*PR  +  XVR*VR  +  XWQ*WQ  +  XW*  W 

1  +  XW*UW  +  XWW*WW  +  XQ*UQ 

1  +  (XDSSDSS*DSDS  +  XDSPDSP*DPDP  +  XDRDR*DRDR  + 

1  +  XDBDB*DBDB)*UU 

1  +  (XWDSP*DSP  +  XWDSS*DSS)*UW  +  (XQDSP*DSP  +  XQDSS*DSS)*UQ 

1  +  (XVDR*UV  +  XRDR*UR)*DR 
C 
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C  HYDRODYNAMIC  LATERAL  FORCE 
C 

DRHS2  =  YPABP*PABP  +  YPQ*PQ 

1  +  YRU*UR  +  YPU*UP  +  YVU*UV  +  YVABVU*UABV  +  YW*W  + 

1  +YWP*WP 

1  +  (YDR  +  YDRN*ETAM1)*UUDR  +  YVABV*VABV  +  YRABR*RABR  + 

1  +  YSTAR*UU 

1  +  YVABVR*VSQRWWW 
C 

C  HYDRODYNAMIC  NORMAL  FORCE 
C 

DRHS3  =  ZQU*UQ  +  ZWU*UW  +  ZABWU*UABW  +  ZVP*VP  +  ZSTAR*UU 

1  +  (ZDSS  +  ZDSSN*ETAM1)  *  UUDS  +  (ZDSP  +  ZDSPN*ETAM1)  *  UUDP 

1  +  ZDB*UUDB 

1  +  ZQABQ*QABQ  +  ZWABW*WABW  +  ZPR*PR  + 

1  +  ZWW*DABS(E(3))*SQRWWW 
C 

C  HYDRODYNAMIC  ROLL  MOMENT 
C 

KBETA  =  BETAS*BETAS*(UU  +  VATS*VATS  +  WATS*WATS) 

DRHS4  =  KPU*UP  +  KRU*UR  +  KQR*QR  +  KWP*WP 

1  +  ((KDR  +  KDRN*ETAM1)*DR  +  KDSS*DSS  +  KDSP*DSP  +  KSTAR)*UU 

1  +  KWR*WR  +  KVU*UV  +  KPABP*PABP  +  KI*VATTAIL  +  ARAPPHI4*K4S 
C      1  +  KBETA*(K4S*DSIN(4.*PHISTERN)  +  K8S*DSIN(8.*PHISTERN)) 
C 

C  HYDRODYNAMIC  PITCH  MOMENT 
C 

DRHS5  =  MQU*UQ  +  MQABQ*QABQ  +  MWABW*WABW 

1  +  MVP*VP  +  MRP*PR  +  MSTAR*UU  +  MDB*UUDB 

1  +  MWU*UW  +  MWABWR*WSQRWWW  +  MABWU*UABW  + 

1  +  MWW*DABS(E(3))*SQRWWW 

1  +  (MDSS  +  MDSSN*ETAM1)*UUDS  +  (MDSP  +  MDSPN*ETAM1)*UUDP 
C 

C  HYDRODYNAMIC  YAW  MOMENT 
C 

DRHS6  =  NPQ*PQ  +  NPU*UP  +  NWP*WP  +  NSTAR*UU 

1  +  NRU*UR  +  NVU*UV  +  NVABVR*VSQRWWW  +  NABVU*UABV  + 

1+  NW*W 

1  +  (NDR  +  NDRN*ETAM1)*UUDR  +  NVABV*VABV  +  NRABR*RABR 
C 

C  CALCULATE  BODY  FORCES,  AND  SUM  FORCES  INDEPENDENT  OF  ALPHA 
OR  BETA 

C  USED  FOR  BOTH  LOW  AND  HIGH  ANGLE  OF  ATTACK  MODELS 
C 
C  BODY  AXIAL  FORCE 
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BRHS1  =  -MASS*(  WQ  -  VR  -  XG*(QQ  +  RR)  +  ZG*PR) 

1  -WTINH20*SINTHETA  +  FXP  +  ARRAY(1) 

1  +  XDRE2*DRE*DRE*URUR 

1  +  (XDSSE2*DSSE*DSSE  +  XDSPE2*DSPE*DSPE)*USUS 
C 

C  BODY  LATERAL  FORCE 
C 

BRHS2  =  -MASS*(UR  -  WP  +  XG*PQ  +  ZG*QR  ) 

1  +  WTINH20*COSTHETA*SINPHI 

1  +  CROSS(2)  +  VORT(2)  +  ARRAY(2) 

1  +  YDRE*(  DELRCR*(-XSTERN*E(6)  -  E(2))  +  DR*URUR) 

1  +YDRE3*DRE3*URUR 
c  1  +  YDRE3*DRE3*URUR 
C 

C  BODY  NORMAL  FORCE 
C 

BRHS3  =  -MASS*(VP  -  UQ  +  XG*PR  -  ZG*(PP  +  QQ)) 

1  +  WTINH20*COSTHETA*COSPHI 

1  +  CROSS(3)  +  VORT(3)  +  ARRAY(3) 

1  +  ZDSSE*(  DELRCS*(XSTERN*E(5)  -  E(3))  +  DSS*USUS) 

1  +  ZDSPE*(  DELRCS*(XSTERN*E(5)  -  E(3))  +  DSP*USUS) 

1  +  (ZDSPE3*DSPE3  +  ZDSSE3*DSPE3)*USUS 
C 

C  BODY  ROLL  MOMENT 
C 

BRHS4  -  MASS*ZG*(UR  -  WP) 

1  -  (ZG*WT  -  ZB*BOUY)*COSTHETA*SINPHI 

1  +  (YG*WT  -  YB*BOUY)*COSTHETA*COSPHI 

1  +  QPROP  +  VORT(4)  +  ARRAY(4) 

1  +  KDRE*(  DELRCR*(XSTERN*E(6)  -  E(2))  +  DR*URUR) 

1  +  KDRE3*DRE3*URUR 

1  +  KDSSE  *  (DELRCS*(XSTERN*E(5)  -  E(3))  +  DSS*USUS) 

1  +  KDSPE  *  (DELRCS*(XSTERN*E(5)  -  E(3))  +  DSP*USUS) 

1  +  (KDSPE3*DSPE3  +  KDSSE3*DSPE3)*USUS 
C 

C  BODY  PITCH  MOMENT 
C 

BRHS5  =  (-LXX  +  IZZ)*PR  +  MASS*(-ZG*(WQ  -  VR)  +  XG*(VP  -  UQ)) 

1  -  (XG*WT  -  XB*BOUY)*COSTHETA*COSPHI 

1  -  (ZG*WT  -  ZB*BOUY)*SINTHETA 

1  +  CROSS(5)  +  VORT(5)  +  ARRAY(5) 

1  +  MDSSE  *  (DELRCS*(XSTERN*E(5)  -  E(3))  +  DSS*USUS) 

1  +  MDSPE  *  (DELRCS*(XSTERN*E(5)  -  E(3))  +  DSP*USUS) 

1  +  (MDSPE3*DSPE3  +  MDSSE3*DSPE3)*USUS 
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C  BODY  YAW  MOMENT 

C 

BRHS6  =  (-IYY  +  D£X)*PQ  -  MASS*(XG*(UR  -  WP)) 
1  +  (XG*WT  -  XB*BOUY)*COSTHETA*SINPHI 
1  +  (YG*WT  -  YB*BOUY)*SINTHETA 

c     1  +  (XG*WT  -  XB*(BOUY))*COSTHETA*SINPHI 

c     1  +  (YG*WT  -  YB*(BOUY))*SINTHETA 
1  +  CROSS(6)  +  VORT(6)  +  ARRAY(6) 
1  +  NDRE*(DELRCR*(-XSTERN*E(6)  -  E(2))  +  DR*URUR) 
1  +  NDRE3*DRE3*URUR 

C 

C  SUM  ALL  FORCES 

C 

RHS(l)  =  BRHS1  +  LOFRAC*DRHSl  +  UPFRAC*ZRHS1 
RHS(2)  =  BRHS2  +  LOFRAC*DRHS2  +  UPFRAC*ZRHS2 
RHS(3)  =  BRHS3  +  LOFRAC*DRHS3  +  UPFRAC*ZRHS3 
RHS(4)  =  BRHS4  +  LOFRAC*DRHS4  +  UPFRAC*ZRHS4 
RHS(5)  -  BRHS5  +  LOFRAC*DRHS5  +  UPFRAC*ZRHS5 
RHS(6)  =  BRHS6  +  LOFRAC*DRHS6  +  UPFRAC*ZRHS6 

C 

RETURN 
END 
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Appendix  B.  Trajectory  input  files. 

This  appendix  presents  the  geometry  limitation  file  and  trajectory  command  files  used  to 
generate  the  simulation  trajectories. 
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Geometry  Limitation  File  GEOM.LIM 

93.57    16.50     0.      XC    YC    ZC    [FT] 

85.275    15.04     0.  XMIN  YMIN  ZMESf  [FT] 

93.576  16.50     0.  XMAX  YMAX  ZMAX  [FT] 

8.00    16.5      0.8125  10.5  D2    HRAD  VRAD  SLEN  [FT] 
0.1     1.625     2.0   1.47     VHOST  [FT/SEC]  DIA  [FT]  RHOC  [SLUG/FTA3] 
ATLMAX 
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Traj_l 

2.0     0.04     VELOCITY  AND  ACCELERATION  EXPONENTS 
XVEL,YVEL,ZVEL  CORRECTION  COMMANDS  [FT/SEC]  [SEC] 
P1A       PIT      P2A      P2T      P3A      P3T 

0.00      0.00     0.00     0.00     0.00     0.00 
0.00     20.00     0.00  600.00     0.00  600.00 
0.00     30.00     0.00  700.00     0.00  700.00 
0.00     35.00     0.00  750.00     0.00  720.00 
0.00     40.00     0.00  800.00     0.00  750.00 
0.00    799.00     0.00  850.00     0.00  800.00 
0.00    899.00     0.00  860.00     0.00  810.00 
0.00    900.00     0.00  870.00     0.00  820.00 
0.00    990.00     0.00  880.00     0.00  850.00 
0.00    999.00     0.00  900.00     0.00  900.00 
ROLL,PITCH,YAW  RATE  COMMANDS  [DEG/SEC]  [SEC] 
P4A       P4T      P5A      P5T      P6A      P6T 
0.00      0.00     0.00     0.00     0.00     0.00 
0.00    600.00     0.00  600.00     0.00    20.00 
0.00    700.00     0.00  700.00     0.00    25.00 
0.00    750.00     0.00  750.00     0.00    28.00 
0.00    800.00     0.00  800.00     0.00  750.00 
0.00    900.00     0.00  850.00     0.00  800.00 
0.00    999.00     0.00  860.00     0.00  810.00 
0.00    999.00     0.00  870.00     0.00  820.00 
0.00    999.00     0.00  880.00     0.00  850.00 
0.00    999.00     0.00  900.00     0.00  900.00 
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Traj_2 

2.0     0.04     VELOCITY  AND  ACCELERATION  EXPONENTS 
XVEL,YVEL,ZVEL  CORRECTION  COMMANDS  [FT/SEC]  [SEC] 

P1A       PIT  P2A      P2T      P3A      P3T 

0.00   0.00  0.00  0.00  0.00  0.00 

0.00  20.00  0.00  600.00  0.00  600.00 

0.00  30.00  0.00  700.00  0.00  700.00 

0.00  35.00  0.00  750.00  0.00  720.00 

0.00  40.00  0.00  800.00  0.00  750.00 

0.00  799.00  0.00  850.00  0.00  800.00 

0.00  899.00  0.00  860.00  0.00  810.00 

0.00  900.00  0.00  870.00  0.00  820.00 

0.00  990.00  0.00  880.00  0.00  850.00 

0.00  999.00  0.00  900.00  0.00  900.00 

ROLL,PITCH,YAW  RATE  COMMANDS  [DEG/SEC]  [SEC] 

P4A   P4T  P5A   P5T   P6A   P6T 

0.00   0.00  0.00  0.00  0.00  0.00 

0.00  600.00  0.00  600.00  0.00  20.00 

0.00  700.00  0.00  700.00  -2.00  25.00 

0.00  750.00  0.00  750.00  0.00  30.00 

0.00  800.00  0.00  800.00  0.00  35.00 

0.00  900.00  0.00  850.00  0.00  36.00 

0.00  999.00  0.00  860.00  0.00  810.00 

0.00  999.00  0.00  870.00  0.00  820.00 

0.00  999.00  0.00  880.00  0.00  850.00 

0.00  999.00  0.00  900.00  0.00  900.00 


91 


Traj_3 

2.0  0.04  VELOCITY  AND  ACCELERATION  EXPONENTS 
XVEL,YVEL,ZVEL  CORRECTION  COMMANDS  [FT/SEC]  [SEC] 

P1A       PIT  P2A      P2T      P3A      P3T 

0.00      0.00  0.00     0.00     0.00     0.00 

0.00  10.00  0.00  600.00  0.00  600.00 

0.00  15.00  0.00  700.00  0.00  700.00 

0.00  20.00  0.00  750.00  0.00  720.00 

-4.00  30.00  0.00  800.00  0.00  750.00 
-4.00  799.00  0.00  850.00  0.00  800.00 
-4.00  899.00  0.00  860.00  0.00  810.00 
-4.00  900.00  0.00  870.00  0.00  820.00 
-4.00  990.00  0.00  880.00  0.00  850.00 
-4.00  999.00  0.00  900.00  0.00  900.00 

ROLL,PITCH,YAW  RATE  COMMANDS  [DEG/SEC]  [SEC] 

P4A       P4T      P5A      P5T      P6A      P6T 

0.00      0.00  0.00     0.00     0.00     0.00 

0.00  600.00  0.00  600.00  0.00  20.00 

0.00  700.00  0.00  700.00  0.00  25.00 

0.00  750.00  0.00  750.00  0.00  28.00 

0.00  800.00  0.00  800.00  0.00  750.00 

0.00  900.00  0.00  850.00  0.00  800.00 

0.00  999.00  0.00  860.00  0.00  810.00 

0.00  999.00  0.00  870.00  0.00  820.00 

0.00  999.00  0.00  880.00  0.00  850.00 

0.00  999.00  0.00  900.00  0.00  900.00 
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Traj_4 

Traj_4  uses  Trajl  with  Vhost  =  0.01  in  Geonxlim 
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Appendix  C.  Model  Fit  Analysis  of  Variance. 
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A/L 

Variance 
Source 

Degrees  of 
Freedom 

Sum  of 
Squares 

Mean  Square 

F  Ratio 

Prob>F 

.68 

Model 

1 

0.04653676 

0.046537 

4.0022 

0.0627 

Error 

16 

0.18604585 

0.011628 

Total 

17 

0.23258261 

.84 

Model 

1 

0.13460215 

0.134602 

2.6603 

0.1339 

Error 

10 

0.50595732 

0.050596 

Total 

11 

0.64055947 

1.0 

Model 

1 

11.940734 

11.9407 

25.7911 

<0.0001 

Error 

70 

32.408566 

0.4630 

Total 

71 

44.349300 

1.16 

Model 

1 

3.2888807 

3.28888 

88.9794 

<0.0001 

Error 

16 

0.5913961 

0.03696 

Total 

17 

3.8802768 

1.32 

Model 

1 

4.7664514 

4.76645 

166.9032 

<0.0001 

Error 

34 

0.9709780 

0.02856 

Total 

35 

5.9374294 

1.46 

Model 

1 

0.46515718 

0.465157 

31.1761 

<0.0002 

Error 

10 

0.14920309 

0.01492 

Total 

11 

0.61436 

Table  C- 1  ANOVA  Table  for  PTF  Effect  on  Launchway  Force 


A/L 

Variance 
Source 

Degrees  of 
Freedom 

Sum  of 
Squares 

Mean  Square 

F  Ratio 

Prob>F 

.68 

Model 

1 

0.04653676 

0.046537 

4.0022 

0.0627 

Error 

16 

0.18604585 

0.011628 

Total 

17 

0.23258261 

.84 

Model 

1 

0.13460215 

0.134602 

2.6603 

0.1339 

Error 

10 

0.50595732 

0.050596 

Total 

11 

0.64055947 

1.0 

Model 

1 

11.940734 

11.9407 

25.7911 

<0.0001 

Error 

70 

32.408566 

0.4630 

Total 

71 

44.349300 

1.16 

Model 

1 

3.2888807 

3.28888 

88.9794 

<0.0001 

Error 

16 

0.5913961 

0.03696 

Total 

17 

3.8802768 

1.32 

Model 

1 

4.7664514 

4.76645 

166.9032 

<0.0001 

Error 

34 

0.9709780 

0.02856 

Total 

35 

5.9374294 

1.46 

Model 

1 

0.46515718 

0.465157 

31.1761 

<0.0002 

Error 

10 

0.14920309 

0.01492 

Total 

11 

0.61436 

Table  C-2  ANOVA  Table  for  PTF  effect  on  Launchway  Moment 
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